{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "                dataset      method  KS statistic       p-value  \\\n",
      "0           Exponential  Fixed Bins      0.096510  8.295925e-08   \n",
      "1           Exponential         GMM      0.042184  7.686409e-02   \n",
      "2           Exponential     Laplace      0.194000  2.265926e-30   \n",
      "3           Exponential  Decimal AR      0.037694  1.478059e-01   \n",
      "4    Square + Student t  Fixed Bins      0.063306  1.323141e-03   \n",
      "..                  ...         ...           ...           ...   \n",
      "115  Square + Student t  Decimal AR      0.067020  5.475213e-04   \n",
      "116     ARIMA Residuals  Fixed Bins      0.130438  9.325848e-12   \n",
      "117     ARIMA Residuals         GMM      0.105188  8.587462e-08   \n",
      "118     ARIMA Residuals     Laplace      0.124500  9.526585e-11   \n",
      "119     ARIMA Residuals  Decimal AR      0.100312  3.987595e-07   \n",
      "\n",
      "     Wasserstein Distance  \n",
      "0                0.115748  \n",
      "1                0.050113  \n",
      "2                0.328035  \n",
      "3                0.046562  \n",
      "4                0.082194  \n",
      "..                    ...  \n",
      "115              0.085739  \n",
      "116              0.279867  \n",
      "117              0.231522  \n",
      "118              0.172670  \n",
      "119              0.237887  \n",
      "\n",
      "[120 rows x 5 columns]\n",
      "Exponential\n",
      "square_student_t_mixture\n",
      "Density\n",
      "ARIMA Residuals\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAFLCAYAAAAteMwoAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/bCgiHAAAACXBIWXMAAA9hAAAPYQGoP6dpAAB1aklEQVR4nO3dd1wUR//A8c8dXQFBUBBEBQsqohRrYu8lGiW2xF5iEkvUxEdsP0uMokl8ErtINBFjIdYYC5bYFWsQNcauQUEp0qVz+/uDhw1HP0Bp8369eMnuzs7Ondx352ZmZxSSJEkIgiAI5YqypAsgCIIgFD8R3AVBEMohEdwFQRDKIRHcBUEQyiER3AVBEMohEdwFQRDKIRHcBUEQyiER3AVBEMohEdyFUkc8VycIRadd0gUQ/jVr1iz27duX63Fzc3MuXLjwFkv09q1btw5dXV3Gjx8PwOrVq1mzZg337t0rcB6zZs3iypUrnDx58k0VUxBKPRHcS5lq1aqxZs2aHI/p6Oi85dK8fStXrmTy5Mny9qBBg2jXrl0JlkgQyiYR3EsZXV1dnJycSroYpYalpSWWlpYlXQxBKHNEm3sZc/v2bRwcHJg1a5a879WrV7Rp04YxY8YgSRJ79+7F3t6egIAABgwYQNOmTenbty++vr5qecXGxuLh4UHXrl1xdHTkvffeY/fu3WppOnfuzKpVq1i+fDnvvPMOTZs2Zdy4cTx9+lQt3bVr1xg+fDjNmjWjZcuWuLu7ExERIR/fu3cvjRs3JiAggCFDhuDo6EinTp3YtGmTnMbe3h6ANWvWyL+vXr1a/h0gLS2NjRs38t5779G0aVOcnJwYOnQoly5dKtobKwjljAjupVBqamqOP5Ik0aRJEz7++GP27duHn58fAPPnz0elUrFs2TIUCoWczyeffEKXLl1Ys2YNtra2TJs2jTNnzgCQmJjIRx99xO+//8748eNZt24drq6uzJ07lw0bNqiVx9vbm8ePH+Ph4cHXX3/N7du3cXd3l49fvXqV0aNHo6+vzw8//MCcOXO4cuUKI0eOJDExUU6nUqmYNm0avXv3ZuPGjbi4uPDNN99w7tw5AHx8fAAYOHCg/HtW3333HevWrWPIkCH8+OOPLF68mKioKKZOnUpCQkIxvPuCUD6IZplSJigoCAcHhxyPzZw5k3HjxjFp0iROnjzJokWLmDBhAidOnGDlypVYWFiopR8xYgSTJk0CoF27dgwYMIC1a9fSoUMH9u7dy/3799m5cyfOzs5ymtTUVNatW8fQoUMxMTEBwNjYmHXr1qGlpQVAYGAgq1evJjIyElNTU1asWIGtrS2enp5ymmbNmtGnTx/27NnDsGHDgPRRMBMnTmTQoEEAuLq6cvz4cU6fPk27du3k5ihLS8tcm6ZCQ0OZPn06I0aMkPfp6ekxZcoU7t27J5q0BOF/RHAvZapVq8b69etzPFajRg0gvWN1+fLlDBo0iLlz5zJgwAB69uyZLf2AAQPk3xUKBd26dWP16tUkJiZy5coVrK2t5cCeoV+/fuzevZuAgAA6dOgAgKOjoxy0AbkNPCEhAX19fQICAhg3bhySJJGamgqAjY0NdevW5cKFC3JwB9Sup6urS9WqVYmPjy/w+7NixQoAIiIiePz4Mf/88w+nTp0CIDk5ucD5CEJ5J4J7KaOrq4ujo2O+6Ro1aoS9vT23b9+mU6dOOaapXr262raZmRmSJBETE0N0dDTVqlXLdo65uTkAMTEx8j4DAwO1NEplemueSqUiJiYGlUqFl5cXXl5e2fLT09NT29bX18+Wlybj2m/dusWiRYu4desWBgYG1KtXDysrK0CMjxeEzERwL6N8fHy4ffs2DRs2ZMmSJbRp0wZjY2O1NFFRUXKwBggPD0dLSwsTExOqVKnCP//8ky3fsLAwAExNTQtUjsqVK6NQKBg9ejR9+vTJdjzrjaEo4uLiGD9+PPb29hw6dAg7OzuUSiVnzpzh6NGjxXYdQSgPRIdqGRQUFMTy5csZOHAgGzZsIDY2liVLlmRLd+LECfl3SZI4duwYrq6u6Orq0qJFC4KCgvD391c758CBA+jo6NC0adMClcXQ0JDGjRvz+PFjHB0d5Z/69euzevVqLl++rNFry/hWkJPHjx8TFRXFyJEjqVevnpz27NmzQPo3CUEQ0omaeymTnJzMjRs3cj3eoEED5s6di4GBATNnzqRKlSpMmzaNpUuX0qNHDzp37iyn/eabb0hKSsLW1pZdu3bx6NEjtmzZAoCbmxvbt29n0qRJfP7559SsWZOTJ0+yZ88eJk+enO1bQF6++OILJkyYwJdffkm/fv1IS0tj8+bNBAQEMHHiRI1ev7GxMX/++SdXr16lefPmasdsbW0xNDRkw4YNaGtro62tzdGjR+Xhm2K0jCD8SwT3UiYsLIwhQ4bkenz69On4+fnxww8/UKVKFSB9VMzvv//O/PnzcXFxkdMuXLgQT09Pnj17RuPGjdm8ebMcMA0MDNi6dSsrVqxg5cqVxMXFYWdnx5IlSxg4cKBGZW7bti2bNm1izZo1fP755+jo6ODg4MBPP/2k8eiVTz/9lHXr1vHxxx9z+PBhtWNGRkasW7eOb775hqlTp1K5cmUaNWrEL7/8wscff8y1a9fUbm6CUJEpJNELVe7s3buX2bNn88cff1CzZs2SLo4gCCVAtLkLgiCUQyK4C4IglEOiWUYQBKEcEjV3QRCEckgEd0EQhHJIBHdBEIRyqFSOc1epVKSmpqJUKtWmsBUEIXeSJKFSqdDW1s7zSd/MyttnrTDvQXlVKoN7amoqt27dKuliCEKZ5OjoiK6uboHSltfPmibvQXlVKoN7xh0361SzgiDkLi0tjVu3bmlUYy1vn7XCvAflVakM7hlfD7W0tMrFH5wgvE2aNK9o+lmLjo6Wp70ozcpDE1NRlerb298XUvjrbAqP/VNLuiiCUOEFBwczcOBAgoODS7ooQgGUypp7hvgYCcQzVoJQKsTFxaFSqYiLiyvpoggFUKpr7oIgCELhiOAuCEI2aRosfKJJWuHtKdXNMt7e3qQkq3BpZY9DezFPtyC8LVpKJe5HtvI4IlTelxoeDYD7kV/QvpreqWpXtTrLe40okTIKeSvVwf3ly5ckJabSulIDHvunkhCb3v5uYKTAzrlUF10QyrzHEaH8Hfpc3taOjscMeBIZQqoqtuQKJhRImYmQCbESr6NF56ogCEJBlIngrlDk3zUgavaC8GZJ2lpI//tXKP3KRAQ0rFwZf39/zp26DoBTiwY4tO+ilkbU7AXhzUqrrEdY98ZIumUibFR4ZeZ/KSY2Rn54onZktRIujSBUTCKwlx1lciikUks8WiwIgpCXMnkbtrAyUWtjN7EQwV4QBCGzMhncFQqlWhu7vqEI7oIgCJmVyeCetYO1Wz9XGtVzLuFSCYIglB5lMriDegdrXFzDEi6NIAhC6VJmg3tmBRkqKQiCUJGUi+AOYqikIAhCZmVyKGR+xBJbgiBUdOWm5p6ZsbExO3bs4OrVqwDUqlWLadOmlWyhBEEQ3qJyWcVVKhVImVZwksRqToIgVDDlMribWxihm1gPXjWAVw0w1xXDJAVByC4pKYk5c+bQvHlz2rZty+bNm3NNe/r0ad5//32cnZ3p27cvf/zxx1ssabqQkBD279+Pu7t7vmnLZbMMQGhwDA/+DgLgnfYuHPH5C/+rf8vHGzrUwW1M85IqniAIpcA333zD7du32bJlC8HBwbi7u2NlZUXPnj3V0t29e5fJkyczc+ZMOnTowPnz55k6dSq7d++mYcM3NxQ7Li6OS5cu4efnh5+fH48fP0ahSG+ZWL58eZ7nltvgntWTBy/wOxsgb1etasrRo0e5c+cOAJaWlnz44YclVTxBEN6y+Ph4du3ahZeXFw4ODjg4OPDgwQO2bduWLbgfPHiQ1q1bM3LkSABq167NyZMnOXLkSLEG9+TkZPz9/bl48SKXLl3i9u3bpKWloVAoMDY2pkOHDri4uODsnH9rRIUJ7lmZWxgRGwTPbhgAYOpqUcIlEgThbbp79y6pqalqgdLV1ZUNGzagUqnURt0NGDCAlJSUbHnExhZ9Rao7d+5w8eJF/Pz8uH79OklJSQDY2NjQt29fXFxccHFxoV69ehrlW2GDO6g33VS3qF7CpRGEkpOWlqa2raWl2YIcWc8vKRnliIuLU3sNurq66OrqqqUNCwvD1NRUbb+5uTlJSUlERUVRtWpVeX/dunXVzn3w4AF+fn4MHTq0yGX+4IMPkCQJJycnRo8ejaOjI05OTpiZmRUp3wod3DMTwyeFiuzWrVvy7wYGBjRu3Fij8+/du0dCQkJxF6vQOnTooFaeyZMnM2XKFLU0CQkJ2QJ+xnZycnKueUdERDBlyhRcXFzo0qXoT8IbGhoSFxfH7du3UalUpKamUqlSJZo3b46Ojk6h8xXBPZOXL1/y4MGDki6GILx1jo6OGtfWM7O3ty/G0hReWloat27d4syZM9lq7lnp6ellC+IZ2/r6+jnmHx4ezpgxY5AkiVWrVhXLA5OXLl3C39+f8+fPc/78eTZv3syPP/6IgYEBrVq1ol27drRr145atWpplK8I7oIgoKWlVaTgXpRz3wRDQ8N8y2RhYUFkZCSpqaloa6eHwrCwMPT19TE2Ns6WPiQkRO5Q9fb2Vmu2KQotLS2aN29O8+bNmTZtGpGRkVy8eJFz585x/vx5zpw5gyRJ1KpVSw70rVq1wsDAIM98RXD/H6VSgYWFBfXr1wfA1ta2hEskCMKb1KhRI7S1tblx4wbNm6cPi75+/TqOjo7ZauTx8fGMHz8epVKJt7c31aq9ufmrTE1N6dOnD3369AHSm7zOnTvHuXPn2L17N9u2bUNHR0etKS0nIrj/j7mFESQawasYABp3bKc2Nt65RSOiVDdFm7wglBMGBgb079+fhQsXsnTpUkJDQ9m8eTMeHh5Aei3eyMgIfX19PD09CQwMZOvWrfIxSG++MTIyeqPltLe3x9raGltbW6ytrTl48GCefQIZRHDPJPPomciw14S+iJHHxpuYmJBUWbTJC0J5Mnv2bBYuXMioUaMwNDRkypQpdO/eHYC2bdvi4eGBm5sbR48eJTExkUGDBqmdP2DAAJYtW1bs5ZIkiVu3bnH+/HkuXLhAQEAAqampaGlp4ejoSNu2bfPNQwT3QiptbYyCIGjOwMCA5cuX5/i057179+TffX1933hZQkJCOHfuHBcuXODixYtER0ejUCgwMzOjb9++tGvXjnfffZcqVaoUKD8R3AupZs2aaot0GxgpsHMWb6cgCJrx8PDg/PnzPHr0CIVCQaVKlXB1daVVq1a88847NGrUqFD5imhUQFk7XO3t7Tl2+Cz+V9Pv7k4tGvCZs1j9SRAEzXh7eyNJEs2aNWPo0KH06dMnx6GbmhLBvYCydrjqJNugkl7Iqz/ZRovpCwRB0NyUKVM4f/48N2/eZPbs2Xz11Vc0bdqUFi1a0KJFC5ycnNDT09M4XxHcNZC1wzUzE5MqahORWVlZMWTIkLdeRkEoS6KjowvchlxeTZw4kYkTJ8ozQF64cIELFy6wdu1aJElCR0cHR0dHWrRoQfPmzXFxcaFy5cr55iuCezHJOhFZvRp2hWqTF+34QkURHBzMqFGj2LJlC1ZWViVdnBJnaGhI165d6dq1KwDPnj3j4sWLnD9/nsuXL+Pv74+npydaWlo0bNiQPXv25JmfiBzFKHPN3tGlPglWEq+jNVsFKiH233MUimIvoiCUGnFxcahUKuLi4kq6KKWSjY0NQ4YMYciQIahUKgICAuRafX4PMIEI7m+MQpF9zonME5MB1KlTh88//zzXPPQqK0RNXig3pDQVCq3852IpaLry5NKlS5w8eRITExOGDx8uT39w9+5dGjZsiFKpxNnZGWdnZyZPnlygG6KIFG+IYeXK+Pv7c+7UdQDadXIlTKX+EFTt2rXZsWMHL1++BKB169YEnEuQR+B06+dKo3rOGtf+BaE0Umgp+WfJQRIDXwEQmJj+7z9LDiLpp09vq1/LjNpz3yuxMpYEX19fpk+fLm8fOnSI/fv3k5KSwoABA+jYsSOrVq1SmyHS0NAw33wr1u3xLYuJjSE4OJjg4GDiXsfKQykzfuzt7UkOt+TZDQOe3TAgPswYlaSSzymOhQAEoTRJDHxFwoNQEh6EonwWi0IC5bNYeV9G4K9INm7cSI0aNfjll19YsGABT5484dChQyiVSoYPH86ZM2fYsmWLxvmKmvtbknUoJaQPp3wd/SLXEThZm3ZyaoMXzTZCaVSQUTBmVGIhnTGk6GO6y7KnT58yadIkXF1dcXV15dChQ1y8eJH+/fszd+5cEhIS8PX1Zfz48RrlKyLBW5S5wxWyB/OssjbtdO/XHH1DZzmYm1go1DpgBaE00GQUjCaBvbwOm9TV1VWbYrhhw4bcvHlT3m7WrBmHDx/WOF/RLFPKZW7akSSJi2f/ZNV/f2TVf3/E78KfJV08QcjmTYyCCQ4OZuDAgfJDg+WJlZWV2uiXevXq8fz5c3k7KSmpQLNAZiWCexmj1o4fl3fNXxDeNPNKRkiqN79+ankeNtmvXz/279+Pv78/kL4q1qtXr3j1Kr3/4ezZs4WaP140y5RhOY3IuRn0SB5u2apVK1QqlTwaB6BGjRrFsqivIAAY6RmgUGpxx+tzXr94CMCrBAkF8PinqUQbpHcUmTXpiJ3bzHzz0zatjEoloVQW7CEPTdKWVnXr1kVLS4sRI0bQsmVLzM3NUSgUzJw5k9jYWG7evJltquGCEMG9jMuoyQPEvW6AsfRv+3uVKlUIeaDHsxvR8r5apnXVpkkQwV4oDq9fPCQu8DYAesD0Wkr0wlRk1LMrWdYtUD5ahnoolQp2+VwmLPTf0WJxryMBBfv3BmBYORCAatWNGDSkVTG+ipIxYcIE+Xc/Pz/594sXLwLg5OSkNlSyoERwL0dymtws82gcgHfauxAbZCpPk2BTxU7t4So7OzsmT5789gsvlCuVtFRFOj8sNJYXwVGZ9iho4TyK2GgFsdFRuZxVNk2aNIn27dtjZGQkr+UKoFAoMDQ0xMTEpFD5iuBezuQ1uVlOaRrY18e4zr+1/Xr16rH3p2vc/espAC6tGqMwfiZq+kKJ09HJe0HosupNVaZEcK/gcqrtvwx6IS8v2KBRLYj9t6aftVmncePG3Lt3j6dPnwKi5i8IhXH58mU8PT2ZP38+derUITExka1bt/Ly5Uv69u2Lk5OTxnmK4C7kW9vPfDxrs04t06qYm5vLwT2/mj+AtbU1gwcPfoOvSBAKJikpiUWLFnHs2DH09fUZO3YsY8eOzTHtnTt3WLBgAffv36devXosWrSIJk2aFLkM165dY9y4cSgUCl69ekWdOnWYNm0ap0+fRqFQ8Ouvv7Jx40batGmjUb4iuAsayxrsdRPrFbjmD9C0rgNHfP7C/+rfAHTs4Ux0RKK8nfWGIG4GwpvyzTffcPv2bbZs2UJwcDDu7u5YWVnRs2dPtXTx8fFMmDCBvn37smzZMnbs2MEnn3zC8ePHqVSpUpHK4OXlRbVq1dixYweWlpY8evSIM2fO0Lt3b+bNm8fEiRPx9PQUwV14+zSp+WekCX0Ro3YDyLqd+Yag6c2gdevWvLirI397AGjoUAe3Mc2L+ZULZVl8fDy7du3Cy8sLBwcHHBwcePDgAdu2bcsW3A8fPoyenh4zZ85EoVAwd+5czp49i6+vL25ubkUqx99//82HH36IpaUlAKdPn0aSJIYPH07VqlV57733WLduncb5iuAulEpZbxia3QyM1b49AFStaqrWV5D1BpDfDaIw3yZ27tzJixcv5G0HBwciIiLkfTl1TmeeJVR8Y3mz7t69S2pqKs7OzvI+V1dXNmzYgEqlQqn89xnPgIAAXF1dUfxvgieFQoGLiws3btwocnCPiopSe0jJz88PIyMjuZ1dT0+vUA9vlcrgLv1vrHbr9o6kpUpY1aqKoVEl2nZ2AijydnHkURbLUJ5fd0xkIhY10j8gBgZ62DWwVLtm0xa10dZRUsMw/QZgpGtDdOUotLTTP6zGJrpo69Qv8PHGde0I8Asm6FkoAPXsrUiIT1XbtrNsCbH/PkZuZWaHiW6yvM/OsiZXr14lMDB93HbDhg2pWdUZRVxoga+Rebu6hTEKw38/PwWRkTY5ORktLS15v5aWFg3MLNFVauV2KgC1qpiRlpZGJetGoJ37Op/61eqQlpaGrq05kk7uD8brWFchLS2N6hZGaOV9aczMjUhLSyMt7d8nZDN+j42NVRtWqKurm23R6bCwMExNTdX2m5ubk5SURFRUFFWrVlVLW69ePfXrm5mpTeFdWObm5jx58gRI/3+4fv06bdq0kW8u58+fp2bNmhrnq5A0+Ut4S5KTkwu00oggCNk5OjpmC2S5Ka+ftUmTJhEZGSlvT548mSlTpqil2b9/PytXruTUqVPyvmfPntG1a1fOnDkjN5MAjBo1CldXV7XFdVauXIm/vz8///xzkcq6aNEifvvtN2bNmsWdO3fYuXMnixcvZtCgQUyfPp0jR47w2WefMXXqVI3yLZU1d21tbRwdHVEqlfLXIEEQ8iZJEiqVSq3Gmp/y9lmTJInExESOHTum1qyS081OT08v24RcGdv6+voFSps1XWF89tlnnDlzhvnz56NQKKhduzZ9+/YF4OXLl/Tu3ZtPP/1U43xLZXBXKpUFrnkIglB45fGzVpBVigAsLCyIjIwkNTVVviGGhYWhr6+vNgVvRtrw8HC1feHh4VSvXr3I5a1evTq//fYbR44cISUlhb59+8o3DW9vb7UVmDQhZoUUBKFCatSoEdra2ty4cUPed/36dfmbTGbNmjXD399f7qOQJIk///yTZs2aFUtZjIyMeP/99/nwww/VbiyFDewggrsgCBWUgYEB/fv3Z+HChdy8eZMTJ06wefNmRo4cCaTX4hMTEwHo2bMnMTExLFmyhIcPH7JkyRISEhLo1atXsZRlw4YNvPPOO9y7l75+8pMnT/j444/p27dvodv0RXAXBKHCmj17Ng4ODowaNYpFixYxZcoUunfvDkDbtm3lFZAMDQ3x9PTk+vXruLm5ERAQwMaNG4v8ABPA3r17WblyJTY2NnLz0KRJk7h48SKSJLFs2TK2bdumcb6lcrSMIAhCRTFkyBAAtm3bhra2NteuXWPEiBFMmjSJyZMnM2PGDJ48ecKePXs0ylfU3AVBEErQP//8Q+/eveVa+7lz55AkSW7yadmypTx3kyZEcBfeurLyZbEg5Swrr0UovRISEjAw+HfupUuXLmFhYUHduukLnKSkpBTq70wE9wK4f/8+06dP591336VJkya0bduWadOmcffu3ZIuWqlx8uRJRo0aRfPmzXF0dKRbt24sWbJEXgcyw65du1i+fHmxXHPWrFl07ty5WPLKLDk5maVLl/L777/nma44X4tQcdWqVYvr19OXyoyIiOCvv/7inXfeASAxMZG9e/dSv359jfMVwT0fDx48YMiQIURFRTFv3jw2b97MzJkzCQ4OZvDgwWrDqCqqffv2MXHiRGxtbfn222/x8vJi1KhRHD16lCFDhhAd/e8yf+vXrycqKqrkClsAoaGhbNmyhdTU1DzTlYXXIpR+7733Hr/99hszZszg008/JTU1lR49egAwaNAg/vrrr0ItkFMqH2IqTX766SdMTU3x8vJSe/Kva9eu9OzZk3Xr1rFx48YSLGHx27t3L7Nnz5aHZeVn7dq19OnTh4ULF8r7WrduTfPmzXn//ffZtWsX48ePf0OlFYSybcyYMQQEBHDo0CEAunfvTseOHYH0icyGDRvGgAEDNM5X1NzzER4eLj/WnVmlSpWYM2dOtnGu27dvp0ePHjRt2pRhw4bh5+eHvb09ly9fBmD16tXY29tnu469vT2rV6+Wt58/f87MmTNp27YtDg4OtGnThpkzZ6rNl9G5c2eWLl3KqFGjaNq0KXPnzgXSZ5mbP38+77zzDo6OjgwePFht4d3ilvEeZdWwYUNmz54tL2jQuXNngoKC2LdvH/b29jx//rzA70d0dDSzZ8+mZcuWtGjRgm+//Tbb/wnAiRMncHNzw9HRkXfffZevv/6a+Ph4+fjq1avp1q0bp0+fpm/fvjRp0oQePXqwf/9+IP1979KlC5A+TC63Zp+cXosgFIauri7r1q3j8OHD/P7776xatUo+tnDhwkIvaylq7vno2LEjZ86cYejQoXzwwQe0bt0aOzs7FApFtjmft27dytdff82IESPo0KED58+fZ9q0aRpfMyEhgZEjR2JqasqCBQswMjLC39+fNWvWoK+vz1dffSWn3bZtG2PGjOHjjz+mcuXKJCUlMWrUKMLDw5k+fTrVq1dnz549jB8/nh9//DHXCf9VKpUcLDP+zdwsoaWllevcIx07duTQoUMkJSXRq1cvWrRogYWFBQCjR4+W061Zs4YJEybQuHFjJk6cWOBHt1UqFePHjycoKAh3d3dMTEz48ccfuXXrlloev//+OzNmzKBv375MmzaNoKAgvv/+ex4+fMhPP/0klz8sLIyvvvqKzz77DGtrazZt2oS7uzuOjo7Y2NiwZs0aJk+ezGeffSaPec6qsK9FEHJja2vLq1evCAgIQKFQUKdOnWzTIGhCBPd8fPTRR4SFhbFp0yY5qJqamtK2bVtGjhxJ06ZNgfQAtH79enr06MG8efMAaNeuHa9fv2bXrl0aXfPp06dYWlqyfPlybGxsgPRmjoCAAK5cuaKW1srKihkzZsjbv/76K3fv3uXXX3+VH41u3749I0aM4Lvvvst1rOycOXPYt2+f2j4HBwf5d29vb1q1apXjuYsXL0alUnHs2DFOnDgBpHcSdenShTFjxsiBvnHjxujq6lK1alWN1oQ8e/YsN2/exMvLi/bt2wPQpk0btVq1JEl89913tGvXju+++07eX6dOHUaPHs2ZM2fkr7oJCQksWbJEvtHVqVOHTp06cebMGcaOHUujRo3k19C4ceMcy1TY1yIIObl//z6LFi2SO1Yhfd6fd999lwULFhRqyl/RLFMAU6dO5dy5c6xYsYKBAwdiaGjI77//zuDBg/H29gbSHxd+9eqV/JU+Q79+/TS+XqNGjdi+fTvW1tY8ffqUM2fOsGnTJh4/fpxtZrqMQJTBz8+PatWq4eDgQGpqKqmpqaSlpdGpUydu376t1rmZ2eTJk9m9eze7d++WF7jO2N69e7daoM/KyMiIVatWceLECebPn0+PHj2IiYnhp59+omfPnvj7+2v8HmR27do1dHR0aNeunbyvUqVKdOjQQd5+/PgxL1++pHPnzvLrTk1NpUWLFhgaGnLhwgW1PDMH5IypXTM33wjC2/L8+XNGjBiBv78/7dq1o1u3bigUCtq3b8+NGzcYOnQoISEhGucrau4FVKVKFd577z3ee+89IH2x3P/85z98++239O3bVx41kXmCf0CutWrqp59+YsOGDURFRWFubk6TJk0wMDAgNjZWLV3Wx5+joqIICwvLNRiHhYVRpUqVbPtr1qwp1w4yFiBwdHTUqMw1a9Zk2LBhDBs2DJVKxYkTJ5g1axaLFy9m7969GuWVWXR0NCYmJtmahTKvXpPx/i9atIhFixZlyyM0NFRtO/O44oxJosSYdaEkrFu3jvj4eLy9vWnevDl+fn4cP34cDw8PXr9+zQcffMCaNWtYvHixRvmK4J6HkJAQPvjgA6ZOncqgQYPUjjVu3Jjp06czadIknj17hqmpKUC2cd1Zh8plBKi0tDR55ZvXr9XXHf39999ZtmwZ//nPf3Bzc5NvGFOnTs13YQUjIyPq1Kmj1jSRWWG+3uXl6NGjLFiwgB07dmBrayvvVyqVdO/enatXr/Lrr7/men5B3g9TU1MiIyPV0oD6e5vRNjlz5kxatmyZ7To53dAEoTQ4f/48Q4YMoXnz7Gv81qxZk6FDh8od/poQzTJ5MDc3R1tbm+3bt5OUlJTt+OPHj9HT06N27drY2tpSo0YNfH191dKcPHlSbTtjrumMdTIBtXa2jG1jY2PGjx8vB/bXr19z/fr1HEeIZNayZUtevHiBmZkZjo6O8s+FCxf48ccf1YJjcahfvz5RUVFs2bIlx+NPnz6lQYMG8nbWqVQL8n60adOG1NRUuT0f0h80ytzUYmdnh5mZGc+fP1d73RYWFqxYsUJe+7QgCvoeZX0tglAYkZGRap+RrKytrbNVGgtC1NzzoKWlxcKFC5k0aRIffPABw4YNo27duiQkJHDhwgW2bdvG1KlT5VrhzJkz+eKLL5g7dy69evXi5s2bbN68WS3PDh064OHhwfz58xk3bhwvXrxg7dq1VK5cWU7TtGlTduzYwbJly+jUqROhoaFs2rSJ8PDwfGugbm5u/PLLL4wZM4ZPP/2UGjVqcPHiRby8vBg+fHiB5od2c3Mr8KK/dnZ2TJgwAU9PT4KDg+nXrx+Wlpa8evWK3377DT8/P3766Sc5vbGxMXfu3OHKlSs0bdq0QO9HmzZtaNu2LfPmzePVq1dYW1vj7e1NREQEZmZmQPr/1fTp05k/fz5aWlp06tSJmJgY1q1bR0hISJ59BlkZGRkB6f0XdevWzXXO7qyvpThW5REqnmrVqvHo0aNcj9+9e7dwo7EkIV+3b9+Wpk+fLrVv315q0qSJ5OLiIg0fPlw6evRotrS+vr5Snz59JAcHB2ngwIHSTz/9JDVo0EC6dOmSnGbfvn1S9+7dJQcHB6lfv37S+fPnpR49ekirVq2SJEmSVCqVtHLlSql9+/aSo6Oj1LVrV2nx4sWSj4+P1KBBA+nhw4eSJElSp06dJHd392xlCA8Pl2bPni21adNGatKkidSjRw/Jy8tLSktLe0PvkCQdPnxYGjNmjNSqVSvJwcFBatOmjTRlyhTp77//Vkv3+++/y+W6evVqgd4PSZKk+Ph46auvvpJatWolOTk5SXPmzJG+/vprqVOnTmr5Hzp0SBowYIDUpEkTqWXLltKnn34q3b17Vz6+atUqqUGDBtnK36BBA7XreXh4SE5OTlKLFi2k5OTkHF9zTq9FEDQ1a9YsqVmzZtLp06clSZKkixcvSvb29tKrV6+kixcvSg4ODtLixYs1zldM+fuGXb58mZEjR+Y5lFAQhIrr2bNnDBkyhIiICNauXUulSpUYO3YsNWvW5NmzZ9SuXRsfHx9MTEw0ylc0GgqCIJQgGxsbtm3bRseOHeVmV4VCQUJCAh999FGhAjuINndBEIQSZ2try4YNG4D0ByJv3LhR5IXLRbOMIAhCOSRq7oIgCCWooGsSZB1WnR8R3AVBEEqQJEnZno5OS0sjOjqalJQUDAwMcHFx0ThfEdwFQRBK0KlTp3Lcr1KpuHbtGjNnzqRFixYa51ukNnd3d3cGDhxYqAvnRaVSkZqailKpzHWaWUEQ1En/W3dAW1u7wE/PlrfPmqbvQVJSEosWLeLYsWPo6+szduxYxo4dm2Pa06dP8/333xMYGEjNmjWZNm1atokC34RffvmFlStXcvXqVY3OK1LN/eDBgxw4cABra2vc3Nzo378/VlZWRckSSJ9HPL85VARByJmjo2OBR1qU189aQd+Db775htu3b7NlyxaCg4Nxd3fHysoq21oNd+/eZfLkycycOVNeq2Hq1Kns3r2bhg0bvqmXAaTPoZTftCM5KVJwP3v2LAcOHGD//v2sWrWKNWvW0Lp1a9zc3OjWrRt6enqFyjfjjuvo6Fjsc6EIQnmVlpbGrVu3NJrzRtPPWkxMTJEWkHjTNHkP4uPj2bVrF15eXjg4OODg4MCDBw/Ytm1btuB+8OBBWrduzciRIwGoXbs2J0+e5MiRI288uPfq1UuemFATRQruZmZmjBkzhjFjxnD37l327t3LoUOHmDFjBkZGRvTu3Rs3N7dc5+bITcbXQy0tLRHcy4EffviBwMBAebtWrVqFWqFKKBhNmlc0+awFBwczatQotmzZUizf0N+kgrwHd+/eJTU1FWdnZ3mfq6srGzZsQKVSqd0gBgwYQEpKSrY8sk7BXVgHDhxg7dq1rFmzhvr16/Pq1StWr17Ny5cvcXNzY9iwYRrnWWwdqg0bNmTOnDm4u7tz9epV9uzZw6+//sqvv/5K3bp1+fDDDxk4cGCha/NC2RUYGCjPES+UTmlpafmmiYmJQaVSERMTU+h1Ct60jNcRFxendrPS1dXN1kwTFhaGqamp2n5zc3OSkpKIiopSW5uhbt26auc+ePAAPz+/Qq9vmtkff/yBu7s7ZmZmJCYmAjBp0iQCAgKoVKkSZ86cYcWKFfTu3VujfIt1tExMTAyHDx/G19eXq1evIkkSTZs2RaFQsHjxYjZt2sSGDRvynN6yPHvsn0pC7L/91wZGCuycy9+Apcw19eLubBfejIK0u2csAn7v3r1sc+6XNh06dCAhIUHenjx5MlOmTFFLk5CQkC3gZ2xnXfEss4iICKZMmYKLi0uxdKhu2bIFOzs7fHx8MDQ05K+//iIgIIBhw4Yxb948JkyYwNatW99+cE9NTeX06dMcOHCA06dPk5ycjKmpKSNGjOCDDz6gfv36APz5559MmDCB2bNn57qOZ3mXECvxOrr8PxCcuaaesQasULoVpM09Yxpme3t7+XNd2mS0uZ85cyZbzT0rPT29bEE8Yzu36ZvDw8MZM2YMkiSxatWqYpnT/9GjR4wdO1Ze2+D06dNIkiRPu925c2e+//57jfMtUnD/6quvOHz4MNHR0SiVStq2bcsHH3xA586d0dZWz9rFxYWOHTvyxx9/FOWSQhmQOaBbWFgQHBxcgqURCiJrm3uaSoVWlsCVEciUSmW+aUuaoaFhvjcrCwsLIiMjSU1NleNVWFgY+vr6OXYah4SEyB2q3t7e2ZbULKzY2Fi16126dImqVavKi7MrFIo8v0nkpkjBffv27dSqVYsxY8bQv3//fNvh6tatS506dYpySaEMMNVqxrNX6WuU6ibWw8oqhLt375ZwqQRNaCmVuB/ZyuOIf9eeTQ1PX1zd/cgvaF9Nn73Qrmp1lvcaUSJlLKpGjRqhra3NjRs35CXurl+/jqOjY7YaeXx8POPHj0epVOLt7a22fm9R1ahRg7///htIX3EtICCArl27Aunj9n19faldu7bG+RYpuG/dujXfNtXExET5K85nn31WlMsJZURwYAQP/g4CoIF9fSzqWKh9jRdNNWXD44hQ/g59Lm9rR8djBjyJDCFVVTyjREqSgYEB/fv3Z+HChSxdupTQ0FA2b96Mh4cHkF6LNzIyQl9fH09PTwIDA9m6dat8DNKbbzJW7iqsbt268fPPP2NkZMSDBw9ITk6Wg/vIkSO5evUq7u7uGudbpOA+a9Ys5s6dm+vENwcOHGDJkiVcvny5KJcRSrnMHag9evRQO2ZuYQSJRvAqRt5nWqdidqiXdZK2FtL//i0vZs+ezcKFCxk1ahSGhoZMmTKF7t27A9C2bVs8PDxwc3Pj6NGjJCYmMmjQILXzBwwYwLJly4pUhk8//RQ/Pz82btwIgJOTk/w5MjIy4uOPP2bMmDEa56tRcI+IiFBb6y8oKIibN2/meOdSqVT88ccfOS4sLZQvmTtQ07/eqndehQbHyDV5gOoWhVgPUihxaZX1COveGEm3/IzwMjAwYPny5SxfvjzbsXv37sm/Z134vjgZGhqya9cuLl68SFpaGm3btpX7C9atW1fofDX6X9LR0WHKlClER6e3vSkUCjw9PfH09MwxvSRJtG/fvtCFE8onY2NjduzYIc+VIR5qKjvKU2AvTTIGpOTk/v37HDt2jMmTJ2uUp0b/U0ZGRixfvhx/f38kScLT05MuXbrkOCxKS0sLMzMz+vXrp1GBhIrh5cuX4sEmQfifR48eER4enuMcMpcuXcLT0xMXFxcUCgU2NjbUrFkz3zw1vg136NCBDh06AOlj14cPH06bNm00zUYQBEEAli9fzk8//ZRnGoVCIc9WWa9ePQ4cOJDvGPsij5YRhKIS8wcJFdmePXto2LAhXbp0yXFOnJyaZYKCgvIddaZRcF+5ciU9e/bE3t5e3s6PQqHg888/1+QyQgVjZWWlNuJGtMELFUlKSgrDhw9n4MCBOR4/ePAgx48ff7Nt7uvXr8fOzk4O7uvXr8/3HBHchYIQk4sJFdX69evlp1Fz0rVr10I92a9RcPf29labHc3b21vjCwr/ymlW0syTi5XXicUEQfhX69atAfDz8+Po0aM8f/4cpVJJnTp1eO+992jatGmhpljWKHK0bNkyz21BM3qVFWrB3MRCUWEmFxME4V/z5s1j9+7dcpu7JEmcO3cOb29vxo0bx3/+8x+N8yxytTA5OZmHDx/KXytu3rzJpk2b0NbWZuTIkRov1FHRZA7m+oZlYw3LzO3jrVq1UuvYqVq1KoHK+JIqmiCUObt372bPnj20a9eOKVOmEBYWxqRJk9i6dSu+vr5s2rSJhg0b0rdvX43yLVJwf/HiBSNGjMDQ0JD9+/cTGhrKqFGj5HmUjx8/zrZt23B0dCzKZSq0zA/71KlTp1T0X2R9IjXzRGE6yTYYG78oyeIJQpny66+/0rhxYzZs2ICWlhZ+fn5A+kSL8+bN48mTJ/zyyy9vN7ivXr2akJAQhg8fDsD+/ftJSEhg+fLlNGvWjPHjx7Nx40ZWr15dlMtUKP7+/pw7dR2Adp1cCVO9VJsbvTSMKslaU79xI0qeXiAyrHQv4iAIpc2DBw+YMmVKrkOCO3XqxIoVKzTOt0jB/eLFiwwbNozRo0cD6Qtmm5ub8/777wMwePDgfAfnVyRnz53lxtX78nb3fs1pWM9JLU1MbIw8/3lMbH2orJ5HaRhV8qZr6mLcu1CRaGtr57iYSIbCzs9VpOD+6tUreeqBhIQEbty4Ic+oBum1uvj4itv+mnW2xMjIWLWFKyRJUqup9x9aNubhyTyl75uoqYtx70JFUqdOHX777TeGDh2abZGjtLQ0Dh06hIODg8b5Fmn5lGrVqvHq1SsgfRhPamoq77zzjnz8wYMHxTqpfVkjSZLa7zk9LpxRUw8ODi7161K+TRnfUB48eCAHeUEojz788ENu377NyJEj1f7Wjx07xtixY7lz5w7jxo3TON8i1dybNm3Kjh07qFWrFl5eXmhra9OpUydSU1M5fvw4u3btolevXkW5RJlWVjoaS9vYeqVSqdauLxb3EMozNzc37ty5w7Zt23j8+DF6enooFAoWLFiArq4uc+bMoVOnThrnW6RP8ZdffsnIkSPlr8wTJ07EzMwMPz8/pk+fjqWlZYVafSlrM0xwYGyZ6Gh822PrlUoFFhb/rs6UdZ3VatWqkRxuKd8YxeIeQnk3b948eanS+Ph4Pv74Y6ytrenYsWO+y5fmpkjB3cbGhgMHDnDx4kVq1KhB06ZNAahfvz7Tpk1j0KBBmJmZFeUSZUp+i1YI6bKuzpTTOquZ2/XNq5mXSDkF4W1q0qSJ/PsXX3xR5PyK/P3byMgo+9Jq5uZ8+umnRc1aKAE5TYnwJmRenale/XrZRgVlllNfRWlrShKEwtq3b1+B0g0YMACA58+fc/XqVXk7N0X+REiSxKVLlwgLC8txonmA/v37F/UywluSdUoEA2MFD0L+4M6dO0D6PBj5zSNd3LKu3NSiRQuaWg8U0zSUYtHR0VSpUqWki1EmzJkzB0mSUCgUaoMwMmTszwjm/v7+zJ49+80G9ydPnjB+/Hi19tLMMgosgvub8abGg2edEiE2yJRnN9Lbv5vWNcbY+O33H2ReucnGxoam1m+9CEIBBQcHM2rUKLZs2VKoCa9yU15vGIsXL9YovZOTE19//XW+6YoU3FesWEFwcDD9+/enSZMmeQ7EFzSXX8fj2xoPnrkJ5U10DGd9nWJ+mrLDvJIRkioNhfLfikZcXBwqlYq4uDi1tFnTaeJN3TBKg9zmcc+NjY1NgUaQFSm4X7p0iQ8//JD58+cXJRshFwXpeCwNT6wWVdbXmXXYaE43uczeVj+BkJ2RngEKpRZ3vD7n9YuHAATFpjfP3vH6nGij9Ca8yjXq0fjjVdnOL2htPLcbhpC7IgX3tLQ0GjVqVFxlEXKgScdjcSmJp2bz+naQNfgbKdUXZM/WTyA6WN+61y8eEhd4GwBVijYK6qAKvU9cZGqu5+RWG1epJJTKgt2xNUlbWjVq1CjHtvasMlfqCqJIn4AmTZrw119/MWjQoKJkU26Ul6lvM89vU1qems0c/B1d6mOX5Zu5mAe/9DDVSWV6rcdU0sp5gEWG3GrjSqWCXT6XCQuNlfdFRr0EYJfPFUxN0pshq1U3YtCQVsVc+rfP1dU1W3BXqVSEh4fz/PlzrKysaNeuncb5Fim4T506lY8//piWLVvSo0ePCj/hU1l5IlUQ3rSsgV3XuBpSmgqFVsFGWoWFxvIiOEreTkxMBhRERSSTGB+V22ll0i+//JLrsYCAAMaMGcO7776rcb5FCu6enp4YGxvz5ZdfMnv2bKpUqZJtmJxCoeDUqVNFuUyZ8aYn1CouWVd/KmsUirc7FFMoOu1Kxii0lPyz5CCJgenzUYWnxKJAQdiyY0g66XOYG7e0pca47E2B+vrGtHAehY6OwVstd0lr1qwZQ4cOxdPTk27duml0bpGC+6NHj9DS0qJGjRpFyUYopMLOwVIWV3/KzLByZbV+gW79XGlUz7mESyUURGLgKxIehALp3UcL6UTlpwkkkL7Aj55N1VzPrWiBPYOtrS3btm3T+LwiBfeTJ08W5XRBQ1lHjdjb2xPyQK9CzsGSuV8gNrbivO7yxlBM0ZGvS5cuYWCg+Y2tWIcUJCcno62t/dafYKwochoy+Dr6hdwUVN2y4k6vLAhl1Q8//JDj/uTkZG7fvs2VK1fkBZA0UeTgHhUVxapVqzh27BgRERFs2rQJHR0dfvzxR9zd3bG1tS3qJYRM8hoyaGFlkm3qADsnMSQwgxguKZRGnp6eeR53dXVlxowZGudbpL/uqKgohgwZwj///IONjY08nCc6OprTp09z8+ZNfHx8xHzcb1HW9vTMAQ1y7kAti6tBFYYYLimURt7e3jnuVyqVVKtWjVq1ahUq3yIF9zVr1hAUFMRPP/2Evb29vApTly5d8PT0ZOrUqaxbtw4PD4+iXEYogqwBLacO1NI4rl0QKooWLVq8kXyL1Dh+8uRJBg8eTJs2bVBkeQa8Q4cODBkyhMuXLxepgIKQn6xDIwsyHYGYskAo74pUcw8NDaVhw4a5Hq9bty7bt28vyiWEN6C8NcPkNDRS39Alz+YoMWWBUN4V6a/ZzMyMoKCgXI/fv38fU1PTolxCeAPKYzNM1qGRBWmOEm3wQnlWpODevn17du7cyaBBg6hcWX1Gq2vXruHj40Pfvn2LVEChaDLXaKF81NTzI55gFYQiBvfJkydz6tQpBgwYgKurKwqFgu3bt7N582YuXLiAkZEREydOLK6yCvnIKahlrtFC+amp5yVrMw2Ip1iFiqdIwd3CwoIdO3bw1VdfcfbsWSRJ4vjx40D62MwFCxZQs2bNYimokL+sQa0i1NJzk/WmJp5iFSqaIvcg2djY4OXlRUxMDE+fPkWSJGrWrImZmVlxlE/QUHlsTy8OoqlGKK3WrFlTqPMmT56c5/FCBXeVSsWBAwf4448/+Ouvv4iIiEChUFC1alWaNGlCz5496dGjh5iGQCg1cvpW8yrktbzt0soeh/adS7KIQgW1du1aeb3pgpIkqfiDe2BgIJ988olcS9fX16dKlSqkpKQQEhJCUFAQx44do169eqxZs4batWtreglBeCOyfqvJvN26UgMxNFIoMYMHD8bFxaVAaf39/fHx8ck3nUZ/vVFRUXz00UdERkYyduxYBg4cqDZ3TFpaGn///Te//fYbO3bsYNSoUfz222/lcsXynIhvKmWXQqEUQyMroKSkJBYtWsSxY8fQ19dn7NixjB07Nse0d+7cYcGCBdy/f5969eqxaNEimjRpUizlaN68Of369StQWi0treIP7ps2bSIiIgIvL68cVwbR0tKiSZMmNGnShG7dujF27Fi8vb2ZMmWKJpcpNfKryWVd9MLY2Pitl1EoHlmbbTp0ac7zhAju3LkDQOPGjbl37x5Pnz4FoE6dOnz++eclVVyhmHzzzTfcvn2bLVu2EBwcjLu7O1ZWVvTs2VMtXXx8PBMmTKBv374sW7aMHTt28Mknn3D8+HEqVapUQqXPm0bB/fTp0/Tq1atASz61bNmS7t27c/r06TIb3DPX5HJqDjt2+Cz+V+8BFXtkSnmRuZlGt5JEbJApz26kz6NtU8WUtLQ0Hjx4kL5tY8MPP/xAYGD6ep7FFezFrJ5vT3x8PLt27cLLywsHBwccHBx48OAB27ZtyxbcDx8+jJ6eHjNnzkShUDB37lzOnj2Lr68vbm5uRSqHQqFApcp7vdnM0tLSCtRKoNFfTlBQEIMHDy5weicnJ86ePavJJUpUXsvPZX1c3cRCgUpSiZEp5Vjm6ZUb2NfHos6/C6XUqFGDkJAQOW316tXZuXMnL16kr5vr4OBARESEvA1gbW2d7+cnr1k9TSwUJMUjgn8xuXv3LqmpqTg7//v8g6urKxs2bEClUqkF0ICAAPlZHkgPyC4uLty4caPIwb1SpUrExMQUOH1sbCz6+vr5ptPoLyMxMREjI6MCpzc2NiY+Pl6TSwDIUwenpaUV+JwzZ84QGRkJgKmpKR06dNDonIYNGxL8SEng0/QPrINTLapVtUYivSy6lSQePwzi4d30mppTy7rYNbCkbWcnAKxqVcXQqFKBtwtzTlG3S+KapaEMxZFH0xa10daxRdko/bmNJvXtsLNMxqLSc3k7PjYZYtO3rczsMNH9dxugcV07AvyCCXqWvsxcPXsrEuJT1bbDY+J4eO9/f2NV6hIXkyL/zWX9m9TRVxDgFyqfX93CGIXhv5+fgshIm5ycrLbAvZaWFg3MLNFV5r3ofa0qZqSlpVHJuhFo6+WaTr9aHdLS0tC1NUfSyb3WqWNdhbS0NKpbGKGV96UxMzciLS1NLU5k/B4bG4u29r/hTVdXF11d9VWfwsLCMDU1Vdtvbm5OUlISUVFRVK1aVS1tvXr11K9vZiZ/kyuK2rVrc+PGDUaOHFmg9Ddu3MDa2jrfdApJg7+Ehg0b8u233xZ4SoEDBw7g7u7O33//XdBLAOl/aLdu3dLoHEEQ0jk6OmYLZLkpr5+1SZMmyRU3SB8TnrV5eP/+/axcuZJTp07J+549e0bXrl05c+YMlpaW8v5Ro0bh6uqq1vS2cuVK/P39+fnnn4tU1tWrV+Pp6cnmzZtp2bJlnmkvXLjAxx9/zLhx4/jyyy/zTFsqv9Npa2vj6OiIUqnUaOynIFRkkiShUqnUaqz5KW+fNUmSSExM5NixY2rNKjnd7PT09EhOTlbbl7Gdtdkjt7QFaR7Jz8iRI9mzZw9jx47lww8/5P3336dRo0Zq36SePHnC7t27+fnnn6latSrjxo3LN1+Ng/u1a9cK3Fzy559/apo9kD6ksKA1D0EQCq88ftYMDQ0LlM7CwoLIyEhSU1PlG2JYWBj6+vrZRr5ZWFgQHh6uti88PJzq1asXubxVqlRhy5YtTJ48mV9++YWtW7eira2NsbExOjo6REdHk5iYiEKhoE6dOnh6emJiYpJvvhoH919//ZVff/21QGk1fepKEAThbWnUqBHa2trcuHGD5s2bA3D9+nX5m0xmzZo1w8vLS45pkiTx559/8umnnxZLWWrXrs1vv/3G0aNH8fX15fbt20RERBAfH0/VqlVxcHCgS5cu9OvXr8DP02jU5v6m5kAQBEEoCfPnz+fPP/9k6dKlhIaG4u7ujoeHB927dycsLAwjIyP09fWJi4ujW7du9OnTh6FDh7Jz5058fX05duxYqR3nrlFwFwRBKE8SEhJYuHAhx44dw9DQkHHjxjF69GgA7O3t8fDwkIc63rx5kwULFvDo0SPs7e1ZtGgRjRs3LsHS500Ed0EQhBI0YsSIAqXbunUrAH5+fqxbt07ezk2pHC0jCIJQUVy/fj3fZxMyHw8PD+fatWv55itq7oIgCCWooFMPZO5IzfoEbY7pi1QqAYAvv/wSe3t7Nm/enO3YrFmzsLe3z/bj7OxM3759+emnn9TSr169Gnt7e3n78uXL8jnnz5/P8foZbYD29vY8f/482/ELFy5gb29f4IfPMsqQ9cfJyYkePXrw/fffk5qaWqC8CirjdV6+fDnPdPb29qxevbpYrw3pX40L+vVYEIqTUqks0E/Wc/IjmmWKKDY2lhMnTtCgQQN8fHwYM2ZMtuGf1apVUxtpJEkS4eHh7Ny5k2XLlqGnp8dHH32U53WUSiW+vr60bds227HDhw/nee6ePXto0KAB9+/f5/r167i6uhbotWWdVjQyMpKDBw+yYcMGUlNT+c9//lOgfArCwcEBHx+fbI94C0JFcODAAdauXcuaNWuoX78+r169YvXq1bx8+RI3Nze6d++ucZ6i5l5EBw8eBGDu3Lk8ffqUS5cuZUujq6uLk5OT/OPs7Ey3bt3YsGEDlpaW7N27N9/ruLi4cPz48RxrzIcPH6ZRo0Y5nhcTE8OJEycYO3Ystra27Ny5s8CvLXOZnZyc6NSpEytWrMDJyalAZdaEoaEhTk5OBX4ARRDKiz/++AN3d3dev35NYmIikD59go+PD1evXmXq1Kn5VuByIoJ7Ee3Zs4c2bdrQunVrateurVHw1NHRwcDAoEAPevXu3ZuoqKhsN4+7d+/y9OlTevXqleN5v//+O6mpqbRr145+/fpx9OhRoqKiClzGnBgaGmYr84kTJ3Bzc8PR0ZF3332Xr7/+Wm3SuMTERBYuXEj79u3lpRg3bdokH8+pWebKlSsMGTKEZs2a0aNHDy5evKh2zdyacrI2sSQmJrJixQq6d+9OkyZNcHFxYcyYMXnOeXThwgUGDx6Ms7MzLVq04LPPPuPRo0eavVGCUABbtmzBzs4OX19fHB0d+euvvwgICGDYsGFcv36ddu3a5TsyJiciuBfBgwcPuHXrFv379wegf//+/PHHH9keUwZITU2Vf5KTk3n+/DkeHh48efJEPj8v9erVo379+vj6+qrtP3ToEC1btqRatWo5nrdnzx7atWuHubk5/fv3JyUlhX379hXo9WUtc0hICF5eXly4cIH3339fTvf7778zadIk7OzsWLt2LZMnT+bAgQNMnDhR7uVfunQpZ8+exd3dnU2bNtGlSxe++eYb9uzZk+O1//rrL8aOHYuRkRGrVq1i5MiRfPHFFwUqd1YzZ85kz549TJgwgc2bNzN79mwePHjAl19+meMohWfPnjFx4kSaNGnC+vXrWbJkCU+ePGHChAkazbstCAXx6NEj3Nzc5G+tp0+fRpIkeXx9586defz4scb5ijb3ItizZw8mJiZ07py+sPKAAQNYvXo1u3fvVnssOSgoCAcHh2zn16lThwULFvDhhx8W6Hq9evXC29ubhQsXynNhHD58ONdHoO/du8dff/3FqlWrALCysqJ169Zy30B+ciqzlZUVU6ZMYcKECUB6/8F3331Hu3bt+O6779Re2+jRozlz5gwdO3bkypUrvPvuu/Tp0weAVq1aUalSJczMzHK8tqenJ2ZmZqxfvx4dHR0gfSrn6dOn51vuzJKTk3n9+jXz5s2jd+/eQPpCMnFxcSxbtozw8PBsN8abN2+SmJjIJ598goWFBQCWlpb88ccfxMfHi6YjoVjFxsaqzWVz6dIlqlatKj8gpVAosk1aVhAiuBdSSkoKBw4coGvXriQmJpKYmEjlypVxdXXl119/ZcKECXKPdrVq1Vi/fj2Q3ga+bt06AgMDWbZsmdpCAfnp3bs3q1at4tKlS7Rt25aAgABCQkLo3r07f/zxR7b0e/bswdjYmObNm8uLAfTo0YMFCxZw6dIlWrdunef1du/eDaQ/xffzzz9z+fJl5s2bR5cuXeQ0jx8/5uXLl3zyySdq/QEtWrTA0NCQCxcu0LFjR1q1asXOnTt5+fIlHTp0oEOHDkyaNCnXa1+/fp1OnTrJgR2ge/fuajPlFYSurq7c/BMSEsKTJ094+vSpPM1rTh+aZs2aoaenx8CBA+nZsyft27enVatWNG3aVKNrC0JB1KhRQ24ifP36NQEBAXTt2hVIrzz5+vpSu3ZtjfMVwb2QTp8+zatXr9i9e7ccBDM7d+6cvGCIrq4ujo6O8jEXFxc++OADPv74Y3bt2qW2yHhebG1tadSokTxq5vDhw7Rt2zbHBcgzbj4xMTG888472Y7v3Lkz3+CeuczNmzdn9OjRTJ06lZ9//lmeaCmj/X7RokUsWrQoWx6hoekLScydOxdLS0sOHDjA4sWLWbx4Mc7OzixcuJCGDRtmOy86OhpTU1O1fdra2tn2FcS5c+dYunQpjx8/pnLlyjRs2FCeDySnZpmaNWvyyy+/sHHjRnbv3o23tzfGxsZ89NFHTJs2TUyGJxSrbt268fPPP2NkZMSDBw9ITk6Wg/vIkSO5evUq7u7uGucrgnsh7dmzBxsbG5YsWaK2X5IkJk+ezM6dO3NdDcrAwIBly5YxZMgQZs+ezY4dOwocMHr37s2mTZtYsGABvr6+zJgxI8d0p06dIjIyksWLF2e76+/YsYMTJ07w6tWrXJtFslIqlXh4eNCnTx9mzZrFoUOH0NPTk79Ozpw5M8eFBjJuPLq6unz22Wd89tlnBAcHc+rUKdatW8eXX37JoUOHsp1nYmKSre9CkiSio6Pl7Yz3LGs7+OvXr6lcuTIAgYGBTJo0ia5du+Lp6YmNjQ0KhYJt27Zx7ty5XF9v06ZNWbNmDcnJyVy/fh0fHx82bNhAw4YNc+28FoTC+PTTT/Hz82Pjxo0A8vMkAEZGRnz88ccFakbNSnSoFkJYWBjnzp2jT58+tGrVSu2ndevW9OzZkzNnzqitsZlV06ZNGTx4MP7+/uzfv7/A1+7VqxdRUVFs2LCB6OhotSaSzPbs2YOlpSWDBg3KVsYRI0aQkpKSa2dmbqytrZk4cSLPnj3Dy8sLADs7O8zMzHj+/DmOjo7yj4WFBStWrODOnTskJibSo0cP+SEvKysrhg0bRp8+feQ1aLNq06YNZ8+eJSEhQd537tw5UlJS5O2Mtu+XL1/K+6Kjo9VGtdy+fZukpCQmTJhArVq15BtCRmDPqeb+888/06lTJ5KTk9HV1aVNmzYsXrwYINfyCkJhGRoaygt1b9iwgW3btsnNjxkVoMIQwb0Q9u/fT2pqqtw5mFX//v1JS0vLd977adOmUaVKFVasWEFcXFyBrm1jY4OjoyOenp506tQpx+lGQ0NDOXfuHL17987xG4Grqyu1atXCx8dHo/U2AUaPHk3t2rXx8vIiKCgILS0tpk+fzs6dO/n666+5cOECR44cYdy4cdy5cwcHBwf09fVxcHBgzZo1bN26lStXruDj48O+ffvkGkpWkyZNIj4+nnHjxnHy5El2797NnDlz1Nrg7e3tqVGjBmvXruX48eOcOHGC8ePHY2BgIKdxcHBAW1ubb7/9lgsXLnDq1CmmTJnC6dOnAXJc47d169aEhYUxadIkzpw5w/nz55k9eza6urp06tRJo/dLEApCqVTStm1bOnTokK1f6f79+4Wabl0E90LYu3cv9evXp0GDBjked3V1pWbNmuzatSvPVatMTU2ZOnUqYWFhrF27tsDX7927NykpKbneXPbv309aWpo8OiQn77//Ps+fP8+zaSInurq6zJkzh8TERJYvXw7AoEGDWLFihbx4wcKFC6lZsyZbt27FxsYGgK+++go3Nzc2b97M2LFjWbduHQMHDmThwoU5XqdOnTr88ssv8s1j3bp1uLu7q/UvaGlpsWrVKszNzfniiy9YsmQJffr0UXuar3bt2qxYsYKQkBA+++wz5s+fD6TPsKdQKHKcgKlhw4Zs2LCBuLg4vvjiCyZPnkxUVBSbN2/Gzs5Oo/dLEAri0aNHXL58GT8/v2w/hw4dYs2aNVy8eBE/P78cpxjJiZg4TBAEoQQtX7482xxTWWWs/gTpz7wcOHAg3/llRIeqIAhCCdqzZw/29vZ07do1x2bU+/fvc+zYMbUV7YKCguRvxbkRwV0QBKEEpaSkMGLECAYOHJjj8YMHD3L8+HGNlysVzTKCIAgl6NKlSzRq1CjH51UgfW6kiIgIrKysNMpXBHdBEIQStnfvXv744w9MTEz49NNP5SaXkydP5jiCpiBEcBcEQShB27dvl5+jgPTpSo4ePUpKSgqtWrWiUaNGeHt7azynkRgKKQiCUIK2b99O/fr1OX36NBs3biQyMpK9e/diYGCAu7s7Dx8+5Mcff9Q431LZoapSqUhNTUWpVIp5PAShgCRJQqVSoa2tXaBl2KD8fdYK8x6UtOfPn/Pll19iYWGBhYUFLVu25Nq1awwbNozRo0cTHBzMyZMnmTZtmkb5lsrgnpqayq1bt0q6GIJQJjk6OqKrq1ugtOX1s6bJe1DSKleurFbWBg0acOXKFbXtXbt2aZxvqQzuGXdcR0fHQnUkCEJFlJaWxq1btzSqsZa3z1ph3oOSZmNjI686BlC/fn21ZSxfv35dqAXpS2Vwz/h6qKWlVS7+4AThbdKkeaWkP2vR0dG5DgEsirLUxDR48GDmzZtH9+7d6dGjB82aNSM6OpoXL15QvXp1jhw5gqWlpcb5lsrgLghF9dg/lYTYfweCGRgpsHMWf+6lSXBwMKNGjWLLli0aj+EuT9LS0tDT02Pq1KnUrl2batWqoVAo+OSTT4iOjiYkJKRQU/6Kv3ahXEqIlXgdLUb5lmZxcXGoVKoCz4haXmVMZqdQKAgMDCQwMBBIX6PZwMCAPn36MGXKFI3zFcFdEAShBHl4eNCuXTuMjIzktZGBIvcbiOAuCIJQgvr37/9G8hXBXRCEt0JKU6HQyr82WtB05cW+ffsKlG7AgAEa5SuCuyAIb4VCS8k/Sw6SGPgKgMDE9H//WXIQST99LV/9WmbUnvteiZWxJMyZMwdJktTmbM+JCO6CIJRaiYGvSHgQCkAS0aCApGcRJJD7imXlXeZ5ZTJIkkRYWBi+vr68fv2aJUuWaJyvCO6CIAglKLd53AEmTJjAgAEDOHv2LK1bt9Yo34rTsCUIwlsTHR2dbxoDdFBI6f8KOdPW1uaDDz7gyJEjGp8rgrsgCMUqODiYgQMHEhwcnGc6MyqxkM6YUektlaxs0tbWJjw8XPPz3kBZBKHEnT13lhtX78vbTi0a4NC+SwmWqOLQ5OEkQ8rG5F4lJSgoiO3bt1OzZk2NzxXBXSiXIiMj1WqOtSOrlWBphOLwpuahKWkdOnTIcX9ycjKRkZEAzJo1S+N8RbOMIAhFIqne/EiXgjb1lEU6Ojpoa2tn+zEyMsLFxYWvv/6a0aNHa5yvqLkLglAkCqUWd7w+5/WLhwAExaoAuOP1OdFG6fVHsyYdsXObWehrlOd5aE6cOPFG8hXBXRCEInv94iFxgbcBiE/SA2oR//IRcZFJAFSyrFuCpauYRHAXBKFY6SvTUCChr9S8uUbbtDIqlYRSWbD52DVJW5rFxcXh5eXF8ePHef78Odra2tjY2NC9e3fGjRuHvr6+xnmK4C4IQoEVpFPTVCeV6bUeU0lLpXH+WoZ6KJUKdvlcJiw0Vt4fGfUSgF0+VzA1SZ8St1p1IwYNaaXxNTJLSkpi0aJFHDt2DH19fcaOHcvYsWNzTHv69Gm+//57AgMDqVmzJtOmTaNLl6KPwIqLi2Po0KE8evQIY2NjGjZsiCRJBAYGsmbNGnx9ffHx8aFSJc2GjIrgLpRJmRfjEAtxvB2aLK5RmMCeWVhoLC+Co+TtuNfpgT48LJbEeM1rsbn55ptvuH37Nlu2bCE4OBh3d3esrKzo2bOnWrq7d+8yefJkZs6cSYcOHTh//jxTp05l9+7dNGzYsEhlWLt2LU+ePGHp0qVq88dIksSuXbtYsGABGzZs4IsvvtAoX/GJEMoksRjH21feOjXj4+PZtWsXXl5eODg44ODgwIMHD9i2bVu24H7w4EFat27NyJEjAahduzYnT57kyJEjRQ7ux44do3///tkmBlMoFAwePJjr169z9OhREdwFQSh/tLX0AMX//i0ed+/eJTU1FWdnZ3mfq6srGzZsQKVSqS2WMWDAAFJSUrLlERsbm22fpsLCwmjatGmux11cXPD19dU4XxHcBUEgLS3/zk+VSiX/mzn921hYW1/fmBbOo9DRMch2LHNZMn6Pi4tTK5euri66uupPw4aFhWFqaqq239zcnKSkJKKioqhataq8v25d9dE+Dx48wM/Pj6FDhxbthQGWlpY8efIk1+PPnj3DzMxM43xFcBcEgVu3bqltZzxYk9mzZ88A+Oeff+RAr6+vj52d3VspY06BHeDevXskJCSo7evQoYPavsmTJ2dbhzQhISFbwM/YTk5OzrUcERERTJkyBRcXl2LpUO3atSvbt2+nadOm9OzZU/7GkJqaysGDB/nll1/o16+fxvmK4C4IAo6Ojuo1cIUCrSxreGYEe1tbWxo0aPA2i5cne3t7+fe0tDRu3brFmTNnstXcs9LT08sWxDO2cxt6GB4ezpgxY5AkiVWrVhV5nVOAzz77jEuXLvHFF19w/Phxvv/+ewD+/vtvZs+eTe3atcUC2UL5lXl0jIlF9nHNO3bs4OrVqwD06NHjrZatPNDS0srWvOJ+ZCuPI0Ll7dTw6P/t/wXtq+nDIdvWacjUd/u8vYLmIKdmIUNDw3ybiywsLIiMjCQ1NVW+cYWFhaGvr4+xsXG29CEhIXKHqre3t1qzTVEYGRmxfft2tm/fjp7ev30Ktra2zJkzhwEDBmBoaKhxviK4C2VC5tEx+obZg/vLly958OABAM2bNwcx22CRPY4I5e/Q5/K2dnQ8ZsCTyBBSVekdibam1UuodEXXqFEjtLW1uXHjxv/+ZuD69es4Ojpmq5HHx8czfvx4lEol3t7eVKtWvBPRZYyxf/XqFQEBASgUCurUqcOIESMKnacI7oIgFIikrYX0v3/LAwMDA/r378/ChQtZunQpoaGhbN68GQ8PDyC9Fm9kZIS+vj6enp4EBgaydetW+RikB2UjI6Mil+X+/fssWrSI69evy/uUSiXvvvsuCxYsKNSUv2JWSKFCKI620YourbIeYd0bk1a5+IYjlrTZs2fj4ODAqFGjWLRoEVOmTKF79+4AtG3blsOHDwNw9OhREhMTGTRoEG3btpV/CrO2aVbPnz9nxIgR+Pv7065dO7p164ZCoaB9+/bcuHGDoUOHEhISonG+ouYulEn+/v6cO5Vey2nXyTXf9MbGxmrt8i1atCAkJITAwPRH2WvVqsW0adPeWHnLC0m3fIUMAwMDli9fzvLly7Mdu3fvnvx7YcaZF9S6deuIj4/H29ub5s2b4+fnx/Hjx/Hw8OD169d88MEHrFmzJseFtPNSvv6nhAojJjZGnts7JrY+VM47vVKpQJL+faJVkiSCgoLkdnpBKCnnz59nyJAhcrt/ZjVr1mTo0KHs379f43xFcBfKhMzL5vUf2l7j880tjCDRCF7FAKCbWA8rqxDu3r1brOUUBE1FRkbmObTU2tqaV69eaZyvCO5CqfTDDz/ITSY9evQgMjJWrqm/fv1aLa1SqcDCwoL69esDULVqVQKV8dnyDA2O4cHfQQDUq18v39q+ILwN1apV49GjR7kev3v3LtWraz4qSQR3oVTIHMxbtWqVrQklrw7RrLVynWQbjI1fvNkCC0IxadWqFT4+Przzzjtq66lKkoSfnx+//vproaY5EMFdKBUCAwPVxqmbajXj2av0x80LEqwz18ojw17nmRay1/Zr165dlOILQqFNnDiRM2fO8Mknn7B27VoqVaqEQqFgyJAhPHv2jNq1azN58mSN8xXBXSiVggMjNArWmspa27du4qI2mkaMnhHeFhsbG7Zt28by5cupUqUKKSkpKBQKEhIS+Oijj/j8888xMTHROF8R3IUiyzo1QFI8ZWIhjaxt8EmVX4rRM0KJsLW1ZcOGDUD6rJs3btzIcT4cTZTOT51QpmSdGiAx7t9tRS7LW+7cuZMXL9KbWho3bvxWpo3VRGkrj1BxKJXKIgd2EMFdeMP0KivUavaQXrvXjrfj2Y0oAGyqmGJlZVWqhiVaWVmpdfKKZhqhrBHBXdBYfjM0ZpV1STx9Q0WZGJaYuZNXEMoaEdwFjeU3Q2NWmacKgMI9hPS2KZVKbGxs5G0xmkYoa0RwF964zFMFQPaHkEqjatWqkRxuKQ/HrNuqudo3ltLcUSwIIIK7UAoU9AnTty3zcMx32rtw8eyf8jcQpxYN+My56EusCcKbIoK7UGSZm136D23Pq5DXatv5KStPmGb+BlI7MvtiDaIDVihNRHAXiixz0Hv9+nW27YLQ9AnTkpbTdAiiA1YoTURwF/K156er3PvrHwA69nDm8Z1XGtXMy6Oc5ocXhNJEBHdBrTmhTp06fP7552rH7/31D35nAwBo0KhWoWrm5U1O88NnfvBJPAQllDQR3AW1RStsbGyyTbcrZJe1n0AvuT7OzpCWlgaAs7Mzq1at4unTp/I5oh1eeJtEcC/nMgdqSG8+aNVwkNqQvsxPhyqVSlQqlZw+v+l2K7LM/QTvtHdBN7Ge2mIgaWkhog1eKDEiuJcxmo61ztx0AKTPNherPvdL5mGI9vb2hDzQ02i6XSFd5mDfwL4+FnX+fV8BtYeiBOFNE8G9jMn6KH9WWZtUMs+LDlCtkmO2oYuZa5w6yTa8jn5RpkauFIfiHmuftdkGwLRO7kupCUJxE8G9GOU59a2xAjun7G+3Wk08lzSZnTt3Dv+r6auyu7ZqiEP7TmrHFZmmYZQkiZfPo+RADenBOmuHaOiLsjUM8U14E2PtM9fkAapbaL5UmiAUlgjuxSivqW/1DXOeHTHrOVlvEDdv3MP/6t9A+jBElaSSA3O36q4c8flL7bi1URv+eaUHiCYVTb3psfZZh0+2atUKlUrFy5cvAbC0tOTDDz8s9usKuUtKSmLRokUcO3YMfX19xo4dy9ixY3NMe+fOHRYsWMD9+/epV68eixYtokmTJm+5xAUngvtblNPsiPk93fnk4Qu1YYhZPXmgflzUwku3ly9fqi0nmBxuybMb0QCYtbAsyaJVSN988w23b99my5YtBAcH4+7ujpWVFT179lRLFx8fz4QJE+jbty/Lli1jx44dfPLJJxw/fpxKlSqVUOnzVm6Cu6bNG0W9RtZml5ymvs0cuLv3a07Dek7Z0hTH051C2ZBTu/6NG/82mzV1tX8rf8dCuvj4eHbt2oWXlxcODg44ODjw4MEDtm3bli24Hz58GD09PWbOnIlCoWDu3LmcPXsWX19f3NzcSugV5K3c/OXk17yRORDntC+nD1LWPDJPHJVTLTvrduZALUlSmZz6Vig++bXrG1aunO1v7N495GY3l1aNURg/486dOwBYW1szePDgt/wqyo+7d++SmpqKs7OzvM/V1ZUNGzagUqnUhgAHBATg6uoq92kpFApcXFy4ceOGCO6ayBi+9/RWEgqFFoYmClKSISk+fX9O2/qGEhLp27qVJB4/DOLh3fRRI04t6xIXkyJv57TPqWVdAvy0CHoWCkA9eyteBMWpHTevUQm7humdYlXM9FAopTy37ZSWtO3sBIBVrarERCbKx3NLY2hU6Y1tv41rlMYylKbXHROZiEWN9EnHDAz0sGuQ+99IFTM9YiIT0dJODyjGJrpo69SnhmH66KfGde0I8AuW/2arWxijMMw+/DUvGWmTk5OzPWHbwMwSXWXeT9rWqmJGWloalawbgbZerun0q9UhLS0NXVtzJJ3cn5vQsa5CWloa1S2MyO8hXzNzI9LS0uQHx+Dfh8hiY2PR1v43vOnq6mZbui4sLAxTU1O1/ebm5iQlJREVFUXVqlXV0tarV0/9+mZmpfo5BoWkyV/CW5KcnMytW7dKuhiCUCY5OjoWeA3O8vpZmzRpEpGRkfL25MmTmTJlilqa/fv3s3LlSk6dOiXve/bsGV27duXMmTNYWv7bBzJq1ChcXV3VpuZYuXIl/v7+/Pzzz2/uhRRBqay5a2tr4+joiFKpVBvaJwhC7iRJQqVSqdVY81PePmuSJJGYmMixY8fUmlVyutnp6emRnJysti9jW19fv0Bps6YrTUplcC+u1b8FQchbefysGRoaFiidhYUFkZGRpKamyjfEsLAw9PX1MTY2zpY2PDxcbV94eDjVq5feZxfEpCGCIFRIjRo1Qltbmxs3bsj7rl+/Ln+TyaxZs2b4+/vLfRSSJPHnn3/SrFmzt1lkjYjgLghChWRgYED//v1ZuHAhN2/e5MSJE2zevJmRI0cC6bX4xMREAHr27ElMTAxLlizh4cOHLFmyhISEBHr16lWSLyFPpbJDVRAE4W1ISEhg4cKFHDt2DENDQ8aNG8fo0aOB9En0PDw85KGON2/eZMGCBTx69Ah7e3sWLVpE48aNS7D0eRPBXRAEoRwqlR2qxSktLY2UlJSSLoYgCG+ZlpYWOjo6JV2MElNug7skSbx8+ZLo6GiNHuoQBKH80NPTw9zcPNvol4qg3Ab36OhooqKiqFatGpUrVy4XY3gFQSgYSZJISUkhOjqaoKD0uXsqWoAvl8FdkiRCQ0MxNjbG3Ny8pIsjCEIJMDAwwMjIiOfPnxMeHl7hgnu5HAqZMd9ERfvPFARBnUKhoEqVKiQlJVW4vrdyGdxTU1MBNHoMWxCE8imjUzXzBGMVQbkM7hlEO7sgCBU1DpTr4C4IglBRVbjgnqZSleprjBgxghEjRhRjaXJ2+fJl7O3tuXz58hu/VkmQ0t78/3NhrzFr1izs7e1z/fH19QXe3t8CQOfOnZk1a1auxzP+XjL/NGnShC5duvDNN9+QkJCQLW15/dsqKypco7SWUon7ka08jgh9I/nbVa3O8l5v5wMp5E6hpeSfJQdJDHz1RvLXr2VG7bnvFfr8atWqsWbNmhyP1alTB4AFCxYUOv83Zf78+Tg4OADpj+7fvXuXVatWERYWxrfffguAg4MDPj4+2Ra3EN6uChfcAR5HhPJ36POSLobwhiUGviLhwZu5iReVrq4uTk5OeaYpjcGxXr16auVu06YNsbGxrF+/ngULFmBoaIihoWG+r0148ypcs0x5sGvXLtzc3HBycqJp06a8//77HDlyRD6+d+9e7O3tCQgIYMCAATRt2pS+ffvKX/dzc+LECT766COcnZ1p0qQJPXv2ZNu2bWppQkNDcXd3p02bNjg7OzN8+HD8/f3l4yqVio0bN9KtWzeaNGlCjx492Lp1a/G+ARVE5mYZb29v7O3t2bt3r3z80qVLNGzYkLVr18r7rl27xvDhw2nWrBktW7bE3d2diIgItXzv3r3LmDFjcHZ2plOnThw4cKBI5cw65Dhrs8zq1avp1q0bp0+fpm/fvvLfxf79+9XO27JlCz179sTR0ZF27dqxcOFC4uLiilS2ikwE9zJm27ZtzJ8/n65du+Lp6cl3332Hrq4uM2bM4OXLl2ppP/nkE7p06cKaNWuwtbVl2rRpnDlzJsd8T58+zaRJk3BwcGDdunWsXr0aGxsbvvrqKwICAgB4/fo1H374IZcvX+Y///kPa9asQU9Pj7Fjx/L06VMAFi5cyKpVq+jXrx8bNmygZ8+eLF26VC0ACelSU1Oz/eQ2VcaIESNo0aIFy5cvJyIigri4OObMmYOTkxOffvopAFevXmX06NHo6+vzww8/MGfOHK5cucLIkSPlqWtDQkIYPnw4sbGxfPvtt0ydOpXvvvuOkJCQApVZpVLJZU1ISODPP//E29ub/v3757lIRlhYGF999RUjR45k48aN1KxZE3d3dx49egTAwYMH+fbbbxk2bBibNm1i0qRJ/PbbbyxevFiTt1TIpEI2y5Rlz549Y9y4cUycOFHeZ21tjZubG9evX6dPnz7y/hEjRjBp0iQA2rVrx4ABA1i7di0dOnTIlu/Dhw8ZMGAAc+fOlfc5OzvTqlUrLl++TLNmzdi3bx9BQUHs27ePRo0aAeDi4kL//v25evUqkiTx66+/8sUXXzBhwgQA2rZti0KhwNPTk48++ghTU9M38r6UNUFBQXLbdWZffvml/N5lplAo8PDwoF+/fnz77bdoaWkRFRXFli1b5IWtV6xYga2tLZ6envK+Zs2a0adPH/bs2cOwYcP4+eefSUtLY+PGjfIC0La2tgwePLhA5c6YDjezmjVrMm3atDzPS0hIYMmSJbRp0wZI71fo1KkTZ86coW7duly5coWaNWsybNgwlEolLVu2pFKlSkRHRxeoXEJ2IriXMRkjGmJiYnj8+DH//POP/PU36xqPAwYMkH9XKBR069aN1atXy7W4zMaPHw+k186fPHlCYGCgvHByRr7Xr1+nZs2acmCH9Ee8jx49CsCOHTuQJInOnTvLD5JB+kiM9evXc/36dbp27Vrk96A8qFatGuvXr8+2P/OizFnZ2NgwY8YMFi9ejCRJeHh4YGNjA6QHz4CAAMaNG4ckSfL7b2NjQ926dblw4QLDhg3j+vXrODk5yYEd0m8AVlZWBSr3okWL5JtScnIyz549Y+PGjQwcOBAfH58888ncDp/xOuPj4wFo3bo1Pj4+uLm50bVrVzp06EDfvn0r7Bj14iCCexkTGBjI/Pnz8fPzQ0dHBzs7Oxo2bAiQ7St91vUdzczMkCSJmJiYbPlGRESwYMECTpw4gUKhoHbt2jRv3lwt36ioKMzMzHItW1RUFIDat4fMCvrVvyLQ1dXF0dFR4/N69+7NsmXLAHj33Xfl/TExMahUKry8vPDy8sp2np6eHpA+oV7NmjWzHa9WrVqBrm9ra6tWbldXV1q2bEnXrl3ZvHkz8+bNy/VcAwMD+feMZewy/rZ69+6NSqVi+/btcrOgtbU1M2bMoHfv3gUqm6BOBPcyRKVSMWHCBHR0dNi9e7e8BuTDhw/57bffsqWPiopSmzgtPDwcLS0tTExMsqWdMWMGjx8/5ueff8bZ2RldXV0SEhL49ddf5TQZkzBl9eeff1KlShW5Y23Lli1Urlw5W7qC1g6F3H399ddUrlwZXV1d5s+fj6enJ4A88+no0aNzvLlmBFZTU9NsCz3DvzfmwrCysqJq1apyv0thvffee7z33nvExsZy/vx5vLy8+M9//oOrqysWFhZFyrsiEh2qZUhkZCRPnjxh4MCBODo6ynPnnD17FkgP/pmdOHFC/l2SJI4dO4arq2uOq91fv36d7t2706pVK/l41nybN2/Os2fPePDggXxeUlISU6ZMYffu3XJNPzIyEkdHR/knIiKClStXFimACHDs2DEOHjzI7NmzmT9/PqdPn2bPnj0AGBoa0rhxYx4/fqz23tevX5/Vq1fLTXetW7fG399f7VvUw4cPefbsWaHL9fz5cyIiIuTx+YUxbdo0uX/IyMiIXr16MXHiRFJTUwkNLZ3DWUs7UXMvhV6+fMnPP/+cbX+DBg2wtrZm27ZtWFpaYmxszLlz5/D29gZQe0oQ4JtvviEpKQlbW1t27drFo0eP2LJlS47XbNq0Kb///jsODg5YWlry559/snHjRhQKhZyvm5sbW7du5bPPPuPzzz/H1NQUb29vUlJS+Oijj7CxsaFfv3783//9H0FBQTRp0oQnT57w/fffU7NmzSJ9+Cu6iIgIFi5cSNu2bXn//fcB6Nq1Kx4eHrz77rtYWlrKHdlffvkl/fr1Iy0tjc2bNxMQECB3wI8aNYrdu3czbtw4pkyZQlpaGt9//32BVyx6+PCh3MQjSRLBwcGsXbsWPT09hg8fXujX17p1axYsWMDy5ctp3749MTExrFmzhjp16sjNjoJmKmRwt6taPf9EJZh3YGAgHh4e2fYPHDiQdevWsWTJEmbNmoWuri716tVj/fr1LF26lGvXrqk9rr5w4UI8PT159uwZjRs3ZvPmzXLtOqtly5axePFieehZnTp1WLRoEQcOHODatWtAeu3wl19+4ZtvvmHx4sWoVCqcnJzw9vaWO/Y8PDzw9PRk586dvHz5EjMzM3r37s20adPkERxvi36t3PsHSnPeOVm0aBEJCQksWrRI3jd//nx69+7N3Llz2bRpE23btmXTpk2sWbOGzz//HB0dHRwcHPjpp5/kzkxTU1N27Ngh/w1VrlyZ8ePHc/jw4QKV46uvvpJ/VyqVmJiY4OTkxLffflukm/fQoUNJSUlh586dbN++HX19fdq0acN//vOfCr1UXlGUywWyExMTefLkCba2tujr66sdS1Op0FK+2daot3GNvOzdu5fZs2fzxx9/5Nh5VhFIaSoUWm/2/+BtXEMourziQXlW4f4y30bQLcnALqR7G0FXBHahNBN/nYIgCOWQCO7lkJubG/fu3auwTTKCIIjgLgiCUC6J4C4IglAOieAuCIJQDongLgiCUA6J4C4IglAOieAuCIJQDongLgiCUA6J4C4IglAOVbjgLqnSysQ1Hj16xOLFi+nRowfNmjXD1dWVoUOHsn379myrHNnb2/Pll1/mmtfgwYOxt7dn9erV8j57e3vs7e3573//m+M5KpWKdu3aZVuUuaxQqd78lEmFvcaIESPk99/e3p6GDRvi7OyMm5sb3t7eav+/xSVj0fSc5uMvbqtXr8be3r5AaYcOHYq9vb28mldWWd+rjPfLxcUFNze3HNcxENJVuFkhFUot7nh9zusXD99I/pVr1KPxx6uKlMfhw4eZPXs2devWZcyYMdja2pKYmMiZM2dYunQp586dY926dfISZEqlklOnTpGUlCRPx5rh+fPn8gLXWSmVSnx9ffniiy+yHbt69WqZnkdbqVSwy+cyYaGxbyT/atWNGDSkVaHPb9y4MQsWLAAgLS2N6Ohozp49i4eHB9euXeOHH36QVysqDh07dsTHxyfb6lwl6fHjx/j7+9OgQQN27txJjx49ckyX+b2C9PcrY1rsmTNnYmJikuO6wBVdhQvuAK9fPCQu8HZJFyNHjx49Yvbs2bRr144ffvhBXpADoEOHDrRq1YrPP/+cI0eOyMuPubi4cO3aNc6ePUu3bt3U8jt8+DCNGjXi77//znatjPPu3LlD48aN1Y4dOnQo1/PKirDQWF4ER5V0MXJkaGiotqYopH8Ls7OzY8mSJRw8eJB+/foV2/WqVq2qtm5qabB3716sra355JNPmDFjBv/88w+1a9fOli6n9wqgffv2tGnThr1794rgnoMK1yxT2v34448olUoWLVqkFtgz9OjRg/79+6vts7GxoUmTJvj6+mZLf/jw4VzXNG3RogXm5ubZzktNTeXYsWO5nie8OcOHD8fCwoKdO3eq7d+1axd9+vShSZMmdOzYkdWrV5OWpt78d+bMGYYOHYqTkxNt27Zl/vz58nq5WZtlZs2axbhx4/Dx8aFr1640bdqUoUOH8uTJE06dOkXfvn1p1qwZgwYNynaD37VrF25ubjg5OdG0aVPef/99jhw5otHrTEtLY//+/XTq1ImuXbtSqVIlfHx8NMpDT08PXV1dsYh2LkRwL2X++OMPWrdunedC1MuXL8+2aHDv3r3lppkMjx8/5u7du7kGaS0tLXr06JEtuPv5+ZGUlETnzp2L8EqEwlAqlbRp04abN2/Kbe+enp783//9H23atGHDhg0MGzYMLy8v/u///k8+79SpU3zyySeYmZnxww8/MGPGDE6cOMH06dNzvZa/vz+//PILs2bNwsPDg0ePHjFhwgQ8PDz45JNP+O9//8uLFy+YMWOGfM62bduYP38+Xbt2xdPTk++++w5dXV1mzJjBy5cvC/w6z549S1hYGP3790dfX59evXqxb98+kpOTs6WVJInU1FT5JykpicePHzN79mxev34tr0wlqKuQzTKlVXR0NNHR0TmuaJO1k02hUKitbNSrVy++/fZbtaaZw4cP4+zsnOfC1L1792bbtm1qTTOHDx+mS5cu2drvhbfD3NyclJQUoqKi0NPTY926dQwZMoR58+YB0LZtW0xMTJg3bx5jxoyR10lt1KgRa9askWuyurq6rFy5MscFsQFev37NDz/8QN26dQG4cuUKO3fu5Oeff6ZNmzYA/PPPPyxfvpyYmBiMjY159uwZ48aNk5ftA7C2tsbNzY3r168X+Nve3r17adCgAY6OjkD6TKa7d+/m6NGj9O3bVy3t1atXcXBwUNunUCho0KABK1eupFOnTgW6ZkUjgnspknWB6wz//PMP3bt3V9tnbW3NyZMn5W0rKyucnJzw9fVVC+7Dhg3L85oZK8v7+vrSuHFjkpOTOXHiBN9++20RX41QWBmLoykUCvz9/UlMTKRz587ZRkkBXLhwARsbG+7cucOUKVPUmih69+6d7RteZlWqVJEDO6TfVACaNWsm7zMxMQGQg/usWbPk7cePH/PPP//Ii2/nVOvOSUREBKdOneLTTz+Vm43q16+PtbU1Pj4+2YK7g4ODvLxgaGgoP/zwAykpKfzwww/Y2dkV6JoVkQjupYipqSmVKlUiKChIbX+NGjXYvXu3vL127Vru37+f7fxevXqxcuVKkpKSePLkCU+fPqVnz555XlOhUNCzZ0951My5c+dQKpW8++67hISEFM8LEzQSEhKCvr4+JiYmREVFATBhwoQc04aGhhIdHY0kSXk25eXE0NAwx/2VKlXK9ZzAwEDmz5+Pn58fOjo62NnZyQtYF3TFzgMHDpCSksLq1avVhucCBAUF8ejRI7WbTuXKleUaPqTffPr168fYsWPZu3dvqesoLi1EcC9lOnfuzKlTp4iLi5M/fLq6ump/3Bm1qax69uzJsmXLOHfuHLdu3cq37T5D79692bJlC3///TeHDx+me/fuYlHiEpKamsrly5dxcXFBS0sLY2NjAL777rscm+vMzc0xNDREoVAQERGhdiwpKYlLly6p1cSLQqVSMWHCBHR0dNi9ezeNGjVCW1ubhw8fajTefM+ePTg7O2frD4iPj2fixIns2LFDboLKibm5OfPnz2fq1KksWbKEFStWFPo1lWeiQ7WUmTBhAqmpqcybNy/Hr7mJiYk8e/Ysx3MtLCxwdXXF19eXI0eOFLj908nJCWtra3777TdOnjwpRsmUIB8fH8LCwvjwww+B9Fqqjo4OISEhODo6yj/a2tr897//5fnz51SuXJlGjRpx6tQptbzOnj3LhAkTiu15hcjISJ48ecLAgQPlMmRcB3JvVszs1q1b3L9/Hzc3N1q1aqX206lTJ1q3bs1vv/1GYmJinvn07NmTdu3acfDgQa5cuVL0F1cOVciae+Ua9Upt3vb29nz77bfMnj0bNzc3Bg4ciL29Pampqfj7+7N7927Cw8MZP358juf36tULDw8PFApFtnb6vPTs2RNvb29MTExo2bJlkV5DaVGtulGpzTsuLo4bN24A6UExMjKS8+fP4+PjQ79+/eT/O1NTU8aPH8/KlSuJi4ujVatWhISEsHLlShQKhdwk8vnnn/PZZ5/xxRdf0L9/f8LDw/nvf/9L165dadCgAbdvF/25DjMzM6ytrdm2bRuWlpYYGxtz7tw5vL29AUhISMg3jz179qCjo5Pr3+b777/PxYsXOXz4MG5ubnnmNWfOHPr168fXX3/Nvn371AYYCBUwuEuqtCI/QVqQayiUhf9D69GjB02aNGHHjh3s3r2boKAgJEnCxsaG3r17M3To0By/okN6kF6yZAkdO3bEyKjgAah3795s2rSJXr16FeuTkSVFpZKK9ARpQa+hVBZujPWdO3cYMmQIkN7vUblyZRo0aMDChQsZNGiQWtpp06ZRrVo1tm/fzo8//kiVKlVo06YNX3zxhfx/3KlTJzZs2MCaNWuYNGkSVatWpW/fvkyZMqVoLzKLdevWsWTJEmbNmoWuri716tVj/fr1LF26lGvXrjFixIhcz01KSuLQoUO8++67uTYtdu/enUWLFrFz5858g7udnR0jRoxg8+bN7Nixg+HDhxflpZU7CqmgvSBlSGJiIk+ePMHW1hZ9ff2SLo4gCCWoosaDsl9FEwRBELIRwV0QBKEcEsFdEAShHBLBXRAEoRwSwV0QBKEcKtfBvRwOBBIEQUMVNQ6Uy+Ce8eh8fHx8CZdEEISS9vr1axQKRYWbUqNcPsSkpaWFiYmJ/Nh1pUqVxIT+glCBZMwBHxMTQ0xMDCYmJhXuCdZy+RATpP/nvnz5Up5VTxCEikdLS4vq1atTpUqVClfBK7fBPUNaWhopKSklXQxBEN4ybW1ttLS0KlxQz1Dug7sgCEJFVC47VAVBECo6EdwFQRDKIRHcBUEQyiER3AVBEMohEdwFQRDKIRHcBUEQyiER3AVBEMqh/wdZdfpG8TIZygAAAABJRU5ErkJggg==",
      "text/plain": [
       "<Figure size 350x300 with 6 Axes>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "import pickle\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "sns.set_style(\"whitegrid\")\n",
    "import matplotlib \n",
    "# matplotlib.rcParams.update({'font.size': 18})\n",
    "\n",
    "ds_map = {\n",
    "    \"Exponential\": \"Exponential\",\n",
    "    \"square_student_t_mixture\": \"Square + Student t\",\n",
    "    \"ARIMA Residuals\": \"ARIMA Residuals\",\n",
    "}\n",
    "\n",
    "#increase vertical spacing\n",
    "fig, ax = plt.subplots(\n",
    "    3, 2, \n",
    "    figsize=(3.5, 3),\n",
    "    gridspec_kw={'width_ratios': [6, 1.5], 'hspace': 0.5, 'wspace': 0.1},\n",
    ")\n",
    "\n",
    "w_dist_map = {\n",
    "    \"Exponent.\": \"Exponential\",\n",
    "    \"Square +\\nStudent-t\": \"Square + Student t\",\n",
    "    \"ARIMA\\nResiduals\": \"ARIMA Residuals\",\n",
    "}\n",
    "\n",
    "csv_fn = \"../precomputed_outputs/uncertainty_quantification/wasserstein.csv\"\n",
    "df = pd.read_csv(csv_fn)\n",
    "df[\"dataset\"] = df[\"dataset\"].apply(lambda x: w_dist_map[x])\n",
    "print(df)\n",
    "\n",
    "samples_fn = \"../precomputed_outputs/uncertainty_quantification/samples.pkl\"\n",
    "with open(samples_fn,'rb') as f:\n",
    "    out_all = pickle.load(f)\n",
    "\n",
    "for i,(ds, datadict) in enumerate(out_all.items()):\n",
    "    print(ds)\n",
    "\n",
    "    samples = datadict['data']\n",
    "    rnn_samples = datadict['RNN-8']\n",
    "    # make seaborn histograms\n",
    "\n",
    "    bins = np.quantile(samples,.005),np.quantile(samples,.99)\n",
    "    ax[i,0].set_xlim(*bins)\n",
    "\n",
    "    sns.histplot(\n",
    "        samples,\n",
    "        ax=ax[i,0],\n",
    "        kde=False,\n",
    "        bins=50,\n",
    "        stat='density',\n",
    "        alpha=.7, \n",
    "        binrange=bins,\n",
    "        label='Ground Truth',\n",
    "        color='k',\n",
    "        edgecolor=None\n",
    "    )\n",
    "\n",
    "    ar_color = '#4f01ecff'# sns.color_palette('Dark2')[3]\n",
    "    sns.histplot(\n",
    "        rnn_samples,\n",
    "        ax=ax[i,0],\n",
    "        kde=False,\n",
    "        bins=50,\n",
    "        stat='density',\n",
    "        alpha=.3, \n",
    "        binrange=bins,\n",
    "        label='Decimal AR',\n",
    "        color=ar_color,\n",
    "        edgecolor=None\n",
    "    )\n",
    "    \n",
    "    if i == 1:\n",
    "        y_label = ax[i,0].yaxis.get_label().get_text()\n",
    "        print(y_label)\n",
    "        ax[i,0].set_ylabel(y_label, labelpad=7, fontsize=14)\n",
    "    else:\n",
    "        ax[i,0].set_ylabel(\"\")\n",
    "\n",
    "    if i!=2:\n",
    "        ax[i,0].set_xticks([])\n",
    "\n",
    "    ax[i,0].set_title(ds_map[ds])#, y=0.65, x=0.27, fontsize=14)\n",
    "    ax[i,0].grid(False)\n",
    "    ax[i,0].set_yticks([])\n",
    "    ax[i,0].set_xticklabels([])\n",
    "\n",
    "    _df = df[df[\"dataset\"] == ds_map[ds]]\n",
    "    sns.barplot(\n",
    "        x='dataset',\n",
    "        y='Wasserstein Distance',\n",
    "        hue='method',\n",
    "        hue_order=[\"Laplace\",\"GMM\",\"Fixed Bins\",\"Decimal AR\"],\n",
    "        errwidth=1,\n",
    "        data=_df,\n",
    "        ax=ax[i,1], \n",
    "        palette=sns.color_palette('Dark2')[:2] + sns.color_palette('Dark2')[3:4] + sns.color_palette('Dark2')[2:3]\n",
    "    )\n",
    "    ax[i,1].get_legend().remove()\n",
    "    ax[i,1].set_title(\"\")\n",
    "    ax[i,1].set_xlabel(\"\")\n",
    "    ax[i,1].set_xticklabels([])\n",
    "\n",
    "    #flip y axis labels to be on the right\n",
    "    ax[i,1].yaxis.tick_right()\n",
    "    ax[i,1].yaxis.set_label_position(\"right\")\n",
    "\n",
    "    if i == 1:\n",
    "        y_label = ax[i,1].yaxis.get_label().get_text()\n",
    "        #flipped horizontally \n",
    "        ax[i,1].set_ylabel(y_label, labelpad=20, fontsize=14, rotation=270)\n",
    "    else:\n",
    "        ax[i,1].set_ylabel(\"\")\n",
    "\n",
    "# plt.subplots_adjust(wspace=0, hspace=0)\n",
    "ax[0,1].legend(\n",
    "    # loc='upper right', \n",
    "    ncol=2, \n",
    "    labelspacing=0.5,\n",
    "    handlelength=1.5, \n",
    "    columnspacing=1.5,\n",
    "    fontsize=12,\n",
    "    bbox_to_anchor=(1.25, -3),\n",
    "    frameon=True,\n",
    ")\n",
    "\n",
    "# plt.tight_layout()\n",
    "plt.subplots_adjust(hspace=0.1)  # adjust space between axes\n",
    "plt.savefig('outputs/histograms.pdf', bbox_inches='tight')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# import pickle\n",
    "# import matplotlib \n",
    "# import matplotlib.pyplot as plt\n",
    "# import seaborn as sns\n",
    "# import numpy as np\n",
    "\n",
    "# sns.set_style(\"whitegrid\")\n",
    "# matplotlib.rcParams.update({'font.size': 18})\n",
    "\n",
    "# with open(\"uncertainty_results.pkl\", 'rb') as f:\n",
    "#     out_alls = pickle.load(f)\n",
    "\n",
    "# samples = out_alls[0][\"square_student_t_mixture\"]['data']\n",
    "\n",
    "# fig, ax = plt.subplots(1, 1, figsize=(5, 3))\n",
    "\n",
    "# bins = np.quantile(samples,.009),np.quantile(samples,.991)\n",
    "\n",
    "# # samples = (samples + 7) / 2\n",
    "# # bins = (bins[0] + 7) / 2, (bins[1] + 7) / 2\n",
    "\n",
    "# sns.histplot(samples,ax=ax,kde=False,bins=100,stat='density',alpha=.7, binrange=bins,label='GT',color='k')\n",
    "\n",
    "# ax.set_ylim(0,0.33)\n",
    "# ax.set_xlabel(\"Value\")\n",
    "\n",
    "# y_ticks = ax.get_yticks()\n",
    "# y_tick_labels = ax.get_yticklabels()\n",
    "\n",
    "# y_ticks = y_ticks[:3]\n",
    "# y_tick_labels = y_tick_labels[:3]\n",
    "\n",
    "# ax.set_yticks(y_ticks)\n",
    "# ax.set_yticklabels(y_tick_labels)\n",
    "# ax.set_xlabel(\"\")\n",
    "\n",
    "# ax.set_xlim(*bins)\n",
    "\n",
    "# ax.spines[['right', 'top']].set_visible(False)\n",
    "# ax.spines['left'].set_bounds(0, 0.12)\n",
    "\n",
    "# plt.grid(False)\n",
    "\n",
    "# color1 = '#a60355'\n",
    "# color2 = '#5400c2'\n",
    "\n",
    "# h1, h2, h3, h4 = np.linspace(0.3, 0.12, num=4, endpoint=True)\n",
    "# x1, x2 = ax.get_xlim()\n",
    "# y1, y2 = ax.get_xlim()\n",
    "# x = (x1 + x2) / 2\n",
    "\n",
    "# plt.text(\n",
    "#     x-14, h1 + 0.05, \n",
    "#     \"0.195 $\\\\rightarrow$ ['1','9','5']\", \n",
    "#     color=\"black\", fontsize=12,\n",
    "#     bbox=dict(\n",
    "#         facecolor='none', edgecolor=color1, boxstyle='round', linewidth=0.5\n",
    "#     )\n",
    "# )\n",
    "# plt.text(\n",
    "#     x+3, h1 + 0.05, \n",
    "#     \"0.537 $\\\\rightarrow$ ['5','3','7']\", \n",
    "#     color=\"black\", fontsize=12,\n",
    "#     bbox=dict(\n",
    "#         facecolor='none', edgecolor=color2, boxstyle='round', linewidth=0.5\n",
    "#     )\n",
    "# )\n",
    "\n",
    "# for i in range(10):\n",
    "#     bin_size = (x2 - x1) / 10\n",
    "#     y = x - 4.5 * bin_size + bin_size * i\n",
    "\n",
    "#     if i == 1:\n",
    "#         plt.plot([x, y], [h1, h2], linewidth=1.5, color=color1)\n",
    "#         plt.text(\n",
    "#             (x+y)/2 - 6, (h1+h2)/2 + 0.01, \n",
    "#             \"P$\\,$('1')\", \n",
    "#             fontsize=11, color='black',\n",
    "#             bbox=dict(\n",
    "#                 facecolor='white', edgecolor=color1, boxstyle='round', linewidth=0.5\n",
    "#             )\n",
    "#         )\n",
    "#     elif i == 5:\n",
    "#         plt.plot([x, y], [h1, h2], linewidth=1.5, color=color2)\n",
    "#         plt.text(\n",
    "#             (x+y)/2 + 1, (h1+h2)/2 + 0.01, \n",
    "#             \"P$\\,$('5')\", \n",
    "#             fontsize=11, color='black',\n",
    "#             bbox=dict(\n",
    "#                 facecolor='white', edgecolor=color2, boxstyle='round', linewidth=0.5\n",
    "#             )\n",
    "#         )\n",
    "#     else:\n",
    "#         plt.plot([x, y], [h1, h2], linewidth=0.25, color='black')\n",
    "    \n",
    "#     for j in range(10):\n",
    "#         bin_size = (x2 - x1) / 100\n",
    "#         z = y - 4.5 * bin_size + bin_size * j\n",
    "#         # plt.scatter(z, 0.15, s=5, linewidth=0.5, edgecolor='black',color='grey')\n",
    "#         if i == 1 and j == 9:\n",
    "#             plt.plot([y, z], [h2, h3], linewidth=1.5, color=color1)\n",
    "#             plt.text(\n",
    "#                 (y+z)/2 + 1.5, (h2+h3)/2 - 0.01, \n",
    "#                 \"P$\\,$('9' | '1')\", \n",
    "#                 fontsize=11, color='black',\n",
    "#                 bbox=dict(\n",
    "#                     facecolor='white', edgecolor=color1, boxstyle='round', linewidth=0.5\n",
    "#                 )\n",
    "#             )\n",
    "#         elif i == 5 and j == 2:\n",
    "#             plt.plot([y, z], [h2, h3], linewidth=1.5, color=color2)\n",
    "#             plt.text(\n",
    "#                 (y+z)/2 + 1, (h2+h3)/2 - 0.01, \n",
    "#                 \"P$\\,$('3' | '5')\", \n",
    "#                 fontsize=11, color='black',\n",
    "#                 bbox=dict(\n",
    "#                     facecolor='white', edgecolor=color2, boxstyle='round', linewidth=0.5\n",
    "#                 )\n",
    "#             )\n",
    "#         else:\n",
    "#             plt.plot([y, z], [h2, h3], linewidth=0.1, alpha=0.5, color='black')\n",
    "\n",
    "#         for k in range(10):\n",
    "#             bin_size = (x2 - x1) / 1000\n",
    "#             zz = z - 4.5 * bin_size + bin_size * k\n",
    "#             # plt.scatter(z, 0.15, s=5, linewidth=0.5, edgecolor='black',color='grey')\n",
    "#             if i == 1 and j == 9 and k == 4:\n",
    "#                 plt.plot([z, zz], [h3, h4], linewidth=1.5, color=color1)\n",
    "#                 # plt.plot([z, zz], [h3, h4], linewidth=0.01, alpha=0.1, color='red')\n",
    "#                 plt.text(\n",
    "#                     (z+zz)/2 - 7, (h3+h4)/2 - 0.01, \n",
    "#                     \"P$\\,$('5' | '19')\", \n",
    "#                     fontsize=11, color='black',\n",
    "#                     bbox=dict(\n",
    "#                         facecolor='white', edgecolor=color1, boxstyle='round', linewidth=0.5\n",
    "#                     )\n",
    "#                 )\n",
    "#             elif i == 5 and j == 2 and k == 9:\n",
    "#                 plt.plot([z, zz], [h3, h4], linewidth=1.5, color=color2)\n",
    "#                 plt.text(\n",
    "#                     (z+zz)/2 + 1, (h3+h4)/2 - 0.01, \n",
    "#                     \"P$\\,$('7' | '53')\", \n",
    "#                     fontsize=11, color='black',\n",
    "#                     bbox=dict(\n",
    "#                         facecolor='white', edgecolor=color2, boxstyle='round', linewidth=0.5\n",
    "#                     )\n",
    "#                 )\n",
    "\n",
    "#         if i == 1 and j == 9:\n",
    "#             plt.scatter(z, h3, s=1, linewidth=0.1, edgecolor=color1,color=color1, zorder=10)\n",
    "#         elif i == 5 and j == 2:\n",
    "#             plt.scatter(z, h3, s=1, linewidth=0.1, edgecolor=color2,color=color2, zorder=10)\n",
    "#         else:\n",
    "#             plt.scatter(z, h3, s=1, linewidth=0.1, edgecolor='black',color='grey', zorder=10)\n",
    "\n",
    "#     if i == 1:\n",
    "#         plt.scatter(y, h2, s=10, linewidth=0.5, edgecolor=color1,color=color1, zorder=10)\n",
    "#     elif i == 5:\n",
    "#         plt.scatter(y, h2, s=10, linewidth=0.5, edgecolor=color2,color=color2, zorder=10)\n",
    "#     else:\n",
    "#         plt.scatter(y, h2, s=10, linewidth=0.5, edgecolor='black',color='grey', zorder=10)\n",
    "\n",
    "# plt.scatter(x, h1, s=20, linewidth=0.5, edgecolor='black',color='grey', zorder=10)\n",
    "\n",
    "# # ycoords = ax.get_label_coords()\n",
    "# # print(ycoords)\n",
    "# ax.set_yticks([])\n",
    "\n",
    "# #set xticks to split xlim into 10 equal parts\n",
    "# ax.set_xticks(np.linspace(*ax.get_xlim(), num=6, endpoint=True))\n",
    "# ax.set_xticklabels([\n",
    "#     '0.0', '0.2', '0.4', '0.6', '0.8', '1.0'\n",
    "# ])\n",
    "\n",
    "# ax.yaxis.set_label_coords(-0.04,0.2)\n",
    "\n",
    "# # ax.xaxis.set_label_coords(0.5,-0.13)\n",
    "# plt.show()\n",
    "\n",
    "# fig.savefig(\"hierarchical.pdf\", dpi=300, bbox_inches='tight')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.3"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
